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This article reviews the basic theoretical aspects of graphene, a one atom thick 
allotrope of carbon, with unusual two-dimensional Dirac-like electronic excitations. 
The Dirac electrons can be controlled by application of external electric and magnetic 
fields, or by altering sample geometry and/or topology. We show that the Dirac 
electrons behave in unusual ways in tunneling, confinement, and integer quantum 
Hall effect. We discuss the electronic properties of graphene stacks and show that 
they vary with stacking order and number of layers. Edge (surface) states in graphene 
are strongly dependent on the edge termination (zigzag or armchair) and affect the 
physical properties of nanoribbons. We also discuss how different types of disorder 
modify the Dirac equation leading to unusual spectroscopic and transport properties. 
The effects of electron-electron and electron-phonon interactions in single layer and 
multilayer graphene are also presented. 
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I. INTRODUCTION 

Carbon is the materia prima for life on the planet and 
the basis of all organic chemistry. Because of the flexi- 
bility of its bonding, carbon-based systems show an un- 
limited number of different structures with an equally 
large variety of physical properties. These physical prop- 
erties are, in great part, the result of the dimensionality 
of these structures. Among systems with only carbon 
atoms, graphene - a two-dimensional (2D) allotrope of 
carbon - plays an important role since it is the basis for 
the understanding of the electronic properties in other 
allotropes. Graphene is made out of carbon atoms ar- 
ranged on a honeycomb structure made out of hexagons 
(see Fig. [I]) , and can be thought as composed o f benzene 
rings stripped out from thei r hydr ogen atoms 1 Pauling! . 
Il972f ). Fullerenes ( Andreonil . l2000l ) are molecules where 
carbon atoms are arranged spherically, and hence, from 
the physical point of view, are zero-dimensional objects 
with discrete energy states. Fullerenes can be obtained 
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from graphene with the introduction of pentagons (that 
create positive curvature defects), and hence, fullerenes 
can be thought as wrapped up graphene . Carb on nan- 
otubes ijCharlier et all . 120071 : ISaito et all . Il998l) are ob- 
tained by rolling graphene along a given direction and 
reconnecting the carbon bonds. Hence, carbon nan- 
otubes have only hexagons and can be thought as one- 
dimensional (ID) objects. Graphite, a three dimen- 
sional (3D) allotrope of carbon, became widely known 
to mankind aft er the invention of the pencil in 1564 
(|Petroskil . [l989h and its usefulness as an instrument for 
writing comes from the fact that graphite is made out 
of stacks of graphene layers that are weakly coupled 
by van der Waals forces. Hence, when one presses a 
pencil against a sheet of paper one is actually produc- 
ing graphene stacks and, somewhere among them, there 
could be individual graphene layers. Although graphene 
is the mother for all these different allotropes and has 
been presumably produced every time someone writes 
with a p encil, it was only isola ted 440 years after its in- 
vention I Novoselov et a/.l l2004l b The reason is that, first, 
no one actually expected graphene to exist in the free 
state and, second, even with the benefit of hindsight, no 
experimental tools existed to search for one-atom-thick- 
flake s among the pencil debris co vering macroscopic ar- 
eas ( Geim and MacDonaldl . [20071 ) . Graphene was even- 
tually spotted due to the subtle optical e ffect it creates on 
top o f a cleverly chosen Si02 substrate I Novoselov et all . 
l2004f ) that allo ws its observation with an ordinar y opti- 
cal microscope (lAbergel et all . 120071 : iBlake et~aH . 120071 : 
ICasiraghi et all . l2007h . Hence, graphene is relatively 
straightforward to make, but not so easy to find. 

The structural flexibility of graphene is reflected in its 
electronic properties. The sp 2 hybridization between one 
s-orbital and two p-orbitals leads to a trigonal planar 
structure with a formation of a cr-bond between carbon 
atoms which are separated by 1.42 A. The c-band is re- 
sponsible for the robustness of the lattice structure in all 
allotropes. Due to the Pauli principle these bands have a 
filled shell and hence, form a deep valence band. The un- 
affected p-orbital, which is perpendicular to the planar 
structure, can bind covalently with neighboring carbon 
atoms leading to the formation of a 7r-band. Since each 
p-orbital has one extra electron, the 7r-band is half-filled. 

Half-filled bands in transition elements have played an 
important role in the physics of strongly correlated sys- 
tems since, due to its strong tight binding character, the 
Coulomb energies are very large, leading to strong col- 
lective effects, magnetism, and insulating beh avior due 
to correlation gaps or Mottness i Phillipsl . l2006h . In fact, 
Linus Pauling proposed in the 1950's that, on the basis 
of the electronic properties of benzene, graphe ne should 
be a resonant valence bond structure (RVB) ijPaulingl . 
Il972h . RVB states have become very popular in the lit- 
erature of transition metal oxides, and particularly i n 
studies of cuprate oxides superconductors ( Mapld . [l998l b 
This point of view, should be contrasted with contem- 
poraneous band structure studies of graphene ( Wallace! . 




Figure 1 (Color online) Graphene (top left) is a honeycomb 
lattice of carbon atoms. Graphite (top right) can be viewed 
a stack of graphene layers. Carbon nanotubes are rolled- 
up cylinders of graphene (bottom left). Fullerenes (C60) are 
molecules consisting of wrapped graphe ne by the introduc- 
tion of pentagons on the hexagonal lattice ijCastro Neto et "all . 
l2006af l. 



Il947h that found it to be a semimetal with unusual lin- 
early dispersing electronic excitations called Dirac elec- 
trons. While most of the current experimental data in 
graphene supports the band structure point of view, the 
role of the electron-electron interactions in graphene is a 
subject of intense research. 

It was P. R. Wallace who in 1946 wrote the first pa- 
pers on the band structure of graphene and sh owed the 
unusual semimetallic behavior in this material (|WallaceL 
Il947h . At that point in time, the thought of a 



purely 2D structure was a mere fantasy and Wallace's 
studies of graphene served him as a starting point 
to study graphite, a very important material for nu- 
clear reactors in the post-World War II era. Dur- 
ing the following years, the study of graphite cul- 
minated with the Slonczewski-Weiss-McClure (SWM) 
band structure of graphite which provided a detailed 
desc ription of the electronic properties in this mate- 
rial (|McCmreL Il957t ISlonczewski and Weissj . Il958h and 



(Bovle and Noziered. 19581: Dillon et aLll977l McClure. 


1958. Il964l Soule et al 


. Il964t ISprv and Schererl. ll96C 


Williamson et all. 1965 


). Interestingly enough, from 



1957 to 1968, the assignment of the electron and hole 
states within the SWM model were the opposite to 
what is accepte d toda y. In 1968, Schroeder et al. 
( Schroeder et all . Il968h established the cur rently ac- 



cepte d location of electron and hole pockets 1 McClure] . 
197lh . The SWM model has been revisited in recent years 
because of its inability to describe the van der Waals-like 
interactions between graphene planes, a problem that re- 
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quires the understanding of many-body effects that go 
beyo nd the band structure description ( Rydberg et all . 
l2003h . These issues, however, do not arise in the context 
of a single graphene crystal but they show up with great 
importance when graphene layers are stacked on top of 
each other, as in the case, for instance, of the bilayer 
graphene. Stacking can change the electronic properties 
considerably and the layering structure can be used in 
order to control the electronic properties. 

One of the most interesting aspects of the graphene 
problem is that its low energy excitations are mass- 
less, chiral, Dirac fermions. In neutral graphene the 
chemical potential crosses exactly the Dirac point. This 
particular dispersion, that is only valid at low ener- 
gies, mimics the physics of quantum electrodynamics 
(QED) for massless fermions except by the fact that 
in graphene the Dirac fermions move with a speed vp 
which is 300 times smaller than the speed of light, 
c. Hence, many of the unusual properties of QED 
can show up in gr a phene b ut at much smaller speeds 
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Castro Neto et al . 2006al IKatsnelson and NovoselovL 



20071 : IKatsnelson et aZ.1 . 12006^ 1 . Dirac fermions behave in 
very unusual ways when compared to ordinary electrons 
if subjected to magnetic fields, lead i ng to new physica l 
pheno mena ijGusvnin and Sharapovl . 120051 : IPeres et all 
l2006d) such as the anomalous integer quantum Hall ef - 
fect (IQHE) measured exp erimentally I Novoselov et all . 
l2005al : IZhang et all l2005f ). Besides being qualitatively 
different from the IQHE observed in S i and GaAlAs 
(heterostructures) devices I Stonel . Il992f ). the IQHE in 
graphene can be observed at room temperature because 
of the large cyclotron e nergies for "relativistic" electrons 
( Novoselov et a/J . l2007l ). In fact, the anomalous IQHE is 
the trademark of Dirac fermion behavior. 

Another particularly interesting feature of Dirac 
fermions is their insensitivity to external electrostatic 
potentials due to the so-called Klein paradox, that is, 
the fact that Dirac fermions can be transmitted with 
probability one through a classi c ally forbidden region 
(ICalogeracos and Dombeyl . Il999t lltzykson and Zubeil 
l2006h . In fact, Dirac fermions behave in a very un- 
usual way in the presence of confining potentials lead- 
ing to the phenomenon of z itterbewequnq. or jitt e ry mo - 
tion of the wavefunction l ltzykson and Zuberl . l2006t l. 
In graphene these electrostatic potentials can be eas- 
ily generated by disorder. Since disorder is unavoid- 
able in any material, there has been great interest in 
trying to understand how disorder affects the physics of 
electrons in graphene and its transport properties. In 
fact, under certain conditions, Dirac fermions are im- 
mune to localization effects observe d in ordinary elec- 
trons ( Lee and Ramakrishnanl . fl985l ) and it has been es- 
tablished experimentally that electrons can propagate 
without scattering over la rge distances of the o rder of 
micrometers in graphene (jNovoselov et all . |2004 ). The 
sources of disorder in graphene are many and can vary 
from ordinary effects commonly found in semiconduc- 
tors, such as ionized impurities in the Si substrate, to 



adatoms and various molecules adsorbed in the graphene 
surface, to more unusual defects such as r ipples associ- 
ated w ith the soft structure of graphene 1 Meyer et all 
l2007al ). In fact, graphene is unique in the sense t hat it 
shares properties of soft membranes ( Nelson et a/1 l20"o~ih 
and at the same time it behaves in a metallic way, so 
that the Dirac fermions propagate on a locally curved 
space. Here, analogie s with problems of quantum grav- 
ity become apparent ( Fauser et all 120071 ) . The softness 
of graphene is related with the fact that it has out-of- 
plane vibrational modes (phonons) that cannot be found 
in 3D solids. These flexural modes, responsible for the 
bending properties of graphene, also account for the lack 
of long range structural order in s oft membranes lead - 
ing the phenomenon of crumpling ( Nelson et all 120041 ). 
Nevertheless, the presence of a substrate or scaffolds that 
hold graphene in place, can stabilize a certain degree of 
order in graphene but leaves behind the so-called ripples 
(which can be viewed as frozen flexural modes). 

It was realized very early on that grap hene should 
also present unusual mes oscopic effects 1 Katsnelsonl . 
l2007al : IPeres et all l2006af ). These effects have their 
origin in the boundary conditions required for the 
wavefunctions in mesoscopic samples with various types 
of ed ges graphene can h ave llAkhmerov and Beenakkei . 



2007: Nakada et al 



199ft IPeres et all 12006c : 
Wakabavashi et all . Il999h . The most studied edges, 



zigzag and armchair, have drastically different electronic 
properties. Zigzag edges can sustain edge (surface) states 
and resonances that are not present in the armchair case. 
Moreover, when coupled to conducting leads, the bound- 
ary conditions for a graphene ribbon strongly affects its 
conductance and the chiral Dirac nature of the fermions 
in graphene can be exploited for applications where one 
can control the valley flavor of the ele ctrons besides 
its c harge, the so-called valleytronics ( Rvcerz et all 
l2007h . Furthermore, when superconducting contacts are 
attached to graphene, they lead to the development of 
supercurrent flow and Andreev process es characteristic 
of su perconducting proximity effect ( Heersche et all . 
l2007h . The fact that Cooper pairs can propagate so well 
in graphene attests for the robust electronic coherence in 
this material. In fact, quantum interference phenomena 
such as weak localizat i on, un iversal conductance fluctu- 
ations ( Morozov et all [20061 ). and the Aharonov-Bohm 
effect in graph ene rings have alr e ady been observe d 



experimentally (jRecher et all 120071 : iRusso et 'all . 120071 ) . 
The ballistic electronic propagation in graphene can be 



used f o r field effect devices such as p-n (ICheianov et al 
2007at [Cheianov and Farkol. l2006t iFqgler et all 12007a 



20071: IWilliams et aZ.l. l2i 
p-n-p (jOssipov et all 



P 



Huard et a/.l.l2007l:lLemme eVa/.l.l2007tlTworzYdlo et al . 



Zhang and Fogleij - boOTt l and 



20071) junctions and as "neutrino " 



billiards ijBerrv and Modragonl . Il987l : iMiao et all 120071) . 
It has also been suggested that Coulomb interactions 
are considerably enhanced in smaller geometries, such 
as gr aphene quantum dots ( Milton Pereira Junior et all . 
2007), leading to unusual Coulomb blockade effects 
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( Geim and NovoselovL l2007h and perhaps to magnetic 
phenomena such as the Kondo effect. The amazing 



transport properties of graphene allow for their use in 



detection 
to spin i 



(Schedin et al. 2007 


• Weh 


ine et al.. 


2007 


ijection (Cho et all 


20071 


Hill et all 


2007 



Because of its unusual structural and electronic flex- 
ibility, graphene can be tailored chemically and/or 
structurally i n many different ways: deposition of 



meta l atoms dCalandra and Mauril . 20071; . . 
2007) o r molecules (ILeenaerts et Schedin et al 

20071 : IWehling et all . l2007f ) on top; intercalation 
(as it is done in graphite inter c alated compounds 
dDresselhaus and Dresselhausl. 120021: iDresselhaus et all . 
Il983l : iTanuma and Kamimural . Il985u ): i ncorporation of 
nitro g en and /and boron i n its structure 1 Martins et all 
120071 : IPeres et all . l2007al ) (in analogy with w hat has 
been done in nanotubes ( Stephan et all \l994 \): using 
different substrates that modify the electronic struc- 



Uchoa et al 



ture ( 


Calizo et al. 2007: Das et al 


. 20071 Faueeras et al. 


2007: Giovannetti et all. 12007: 


Varchon et all 2007 


Zhou et al.. 200711. The control of eraphene properties 



can be extended in new directions allowing for creation 
of graphene-based s ystems with magnetic and supe rcon- 
ducting properties (jUchoa and Castro Netd . l2007h that 
are unique in their 2D properties. Although the graphene 
field is still in its infancy, the scientific and technolog- 
ical possibilities of this new material seem to be un- 
limited. The understanding and control of the proper- 
ties of this material can open doors for a new frontier 
in electronics. As the current status of the experiment 
and potential applica t ions h ave recently been reviewed 
( Geim and NovoselovL 120071 ). in this article we mostly 
concentrate on the theory and more technical aspects of 
electronic properties of this exciting new material. 



II. ELEMENTARY ELECTRONIC PROPERTIES OF 
GRAPHENE 

A. Single layer: tight-binding approach 

Graphene is made out of carbon atoms arranged in 
hexagonal structure as shown in Fig. [2 The structure 
is not a Bravais lattice but can be seen as a triangular 
lattice with a basis of two atoms per unit cell. The lattice 
vectors can be written as: 



ai = ^(3,V3) 



a 2 = -(3,-V3), (1) 



where a ~ 1.42 A is the carbon-carbon distance. The 
reciprocal lattice vectors are given by: 



2tt 
3a 



(1,V3), 



2tt 
3a 



(1,-V3). (2) 



Of particular importance for the physics of graphene are 
the two points K and K' at the corners of the graphene 



Brillouin zone (BZ). These are named Dirac points for 
reasons that will become clear later. Their positions in 
momentum space are given by: 



K ^ 



2?r 2tt 



3a ' 3\/3c 



K'= [tL 



2tt 

3a ' 



2tt 



(3) 



The three nearest neighbors vectors in real space are 
given by: 



«5 1 = |(1,V3) 



S 3 = -a(l,0) 
(4) 

while the six second-nearest neighbors are located at: 
S[ = ±ai, <$2 = ± a 2, £3 = ±(02 - ai). 

The tight-binding Hamiltonian for electrons in 
graphene considering that electrons can hop both to near- 
est and next nearest neighbor atoms has the form (we use 
units such that h = 1): 



- * E ( a U«*J + b U b "d + h - c -) ' ( 5 ) 



where a, lCT (aj CT ) annihilates (creates) an electron with 
spin a (a =T,i) on site Ri on sublattice A (an equiva- 
lent definition is used for sublattice B), t (rts 2.8 eV) is the 
nearest neighbor hopping energy (hopping between dif- 
ferent sublattices), t' x is the next nearest neighbor hop- 
ping energy (hopping in the same sublattice). The en- 
ergy bands der ived from this Hamiltonian have the form 
(|Wallacd . fl947h : 



£ ± (k) = ±V3 + /(k)-i'/(k). 



/(k) = 2 cos (v^fc 



- 4 cos k v a cos —k r a 

\ 2 V J \2 , 

(6) 

where the plus sign applies to the upper (it) and the 
minus sign the lower (n*) band. It is clear from (O that 
the spectrum is symmetric around zero energy if t' = 
0. For finite values of t' the electron-hole symmetry is 
broken and the ir and ir* bands become asymmetric. In 
Fig. [3] we show the full band structure of graphene with 
both t and t'. In the same figure we also show a zoom 
in of the band structure close to one of the Dirac points 
(at the K or K' point in the BZ). This dispersion can 
be obtained by expanding the full band structure, eq.((6|), 



1 The value o f t' is not well known but ab initio calculations 
IjReich et all 120021 1 find 0.02t < t' < 0.2t depending on the 
tight-binding parameterization. These calculations also include 
the effect of a third nearest neighbors hopping, which has a value 
of around 0.07 eV. A tight bind ing fit to cyclotron resonance ex- 
periments jDeacon et aij . l2007l ) finds f R3 0.1 eV. 



Figure 2 (Color online) Left: Lattice structure of graphene, 
made out of two interpenetrating triangular lattices (ai and 
a,2 are the lattice unit vectors, and Si, i = 1, 2, 3 are the near- 
est neighbor vectors); Right: corresponding Brillouin zone. 
The Dirac cones are located at the K and K' points. 

close to the K (or K') vec tor, eq.([3]), as: k = K + q, with 
|q| < |K| (|Wallacd . fl947l l: 

E ± (q)^±v F \q\+0((q/K) 2 ), (7) 

where q is the momentum measured relatively to the 
Dirac points and vf represents the Fermi velocity, given 
by vf = Ma/2, with a value vf — 1 x 10 6 m/s. This 
result was first obtained by Wallace ( Wallace! . Il947h . 

The most striking difference between this result and 
the usual case, e(q) = q 2 /(2m) where m is the electron 
mass, is that the Fermi velocity in ([7]) does not depend 
on the energy or mo mentum: in the usual case we have 
V = k/m = y / 2E/m and hence the velocity changes sub- 
stantially with energy. The expansion of the spectrum 
around the Dirac point including t' up to second order in 
q/K is given by: 

E ± ( q ) c 3t'±vM\-( ^ ± ^ sin(30 q )) |q| 2 , (8) 



where 



arctan | — 
q y 



(9) 



is the angle in momentum space. Hence, the presence of t' 
shifts in energy the position of the Dirac point and breaks 
electron-hole symmetry. Notice that up to order (q/K) 2 
the dispersion depends on the direction in momentum 
space and has a three fold symmetry. This is the so-called 
trigonal warping of the electronic spectru m l Ando et all 
ll998l : lDresselhaus and Dresselhausl . l2002h . 



1. Cyclotron mass 

The energy dispersion |(7j) resembles the energy of 
ultra-relativistic particles; these particles are quantum 
mechanically described by the massless Dirac equation 
(see section Hl.Bl for more on this analogy). An immedi- 
ate consequence of this massless Dirac-like dispersion is 
a cyclotron mass that depends on the electronic density 



Figure 3 (Color online) Left: Energy spectrum (in units of 
t) for finite values of t and t' , with t =2.7 eV and t' = 0.2t. 
Right: zoom- in of the energy bands close to one of the Dirac 
points. 



as its square root ( Novoselov et all . l2005aUZhang et all 
l2005h . The cyclotron mass is defined, within the se mi- 
classical approximation (jAshcroft and Merminl . Il97fj ) , as 



2tt 



dA(E) 
dE 



(10) 



E=E F 



with A(E) the area in k— space enclosed by the orbit and 
given by: 



A{E) 

Using l[TTj) in iflOj) one obtains 



nq(E) 2 = tt^- 



E F 



(11) 



(12) 



The electronic density, n, is related to the Fermi momen- 
tum, kp, as kp/n = n (with contributions from the two 
Dirac points K and K ' and spin included) which leads 
to: 



vf 



(13) 



Fitting (fl~3|) to the experimental data (see Fig[4]) pro- 
vides an estimation for the Fermi velocity and the 
hopping parameter as vf « 10 6 ms^ 1 and t w 3eV, 
respectively. The experimental observation of the 
\fn dependence of the cyclotron mass provides evi- 
dence for the e xistence of massles s Dirac quasiparti- 



cles in graphene ^Deacon et al\. 12 007 



iNovoselov et all . l2005at IZhang et al 



I Jiang et a/.l . l2007at 



2005? ) - the usual 



parabolic (Schrodinger) dispersion implies a constant cy- 
clotron mass. 



2. Density of states 

The density of states per unit cell, derived from ((6|), is 
given in Fig. [5] for both t' = and t' ^ 0, showing in 




n (10 12 crrr 2 ) 



Figure 4 (Color online) Cyclotron mass of charge carriers in 
graphene as a function of their concentration n. Positive and 
negative n correspond to electrons and holes, respectively. 
Symbols are the experimental data extracted from tempera- 
ture dependence of the SdH oscillations; solid curves is the 
best fit to Eq. lfl3ll. mo is t he free electron mass. Adapted 
from iNovoselov et ali l2005al . 



both cases a semime tallic behavior i Bena and Kivelsorl 
l2005l : IWallacelll947l ). For t 1 = it is possible to derive 
an analytical expression for the density of states per uni t 
cell, which has the form ( Hobson and Nierenbergl . ll953f l: 




-t<E<t 



-3i < E < -tV t < E < 3i 



-t < E < t 



(f) "I 



-M<E<-tVt < E<3t 
(14) 

where F(tt/2,x) is the complete elliptic integral of the 
first kind. Close to the Dirac point the dispersion is ap- 
proximated by ([7]) and the expression for the density of 
states per unit cell is given by (with a degeneracy of four 
included): 



where A c is the unit cell area given by A c = 3\/3a 2 /2. It 
is worth noting that the density of states for graphene is 
very dif ferent from t he den sity of states of carbon nan- 
otubes ijSaito et all . Il992al lbl). The latter show 1/VE 
singularities due to the ID nature of their electronic 
spectrum, which comes about due to the quantization 
of the momentum in the direction perpendicular to the 
tube axis. From this perspective, graphene nanoribbons, 
which also have momentum quantization perpendicular 
to the ribbon length, have properties very similar to car- 
bon nanotubes. 




Figure 5 (Color online) Density of states per unit cell as a 
function of energy (in units of t) computed from the energy 
dispersion (J6}, t' — 0.2t (top) and for t' = (bottom). Also 
shown is a zoom in of the density of states close to the neu- 
trality point of one electron per site. For the case t' = the 
electron-hole nature of the spectrum is apparent and the den- 
sity of states close to the neutrality point can be approximated 
by p(e) oc |e|. 



B. Dirac fermions 

Let us consider Hamiltonian with t' = and con- 
sider the Fourier transform of the electron operators: 



1 



(16) 



where N c is the number of unit cells. Using this transfor- 
mation, let us write the field a n as a sum of two terms, 
coming from expanding the Fourier sum around K' and 
K. This produces an approximation for the representa- 
tion of the field a n as a sum of two new fields, written 
as 



a„ 
b„ 



-iKR 



-iK' R n 



h, n + e 



-iK'R r 



b2,n ■ 



(17) 



P(E) 



2A r \E\ 



7T Vf 



(15) 



where the index i = 1 (i = 2) refers to the K (K') 
point. These new fields, ai_ n and bi_ n are assumed to 



7 



vary slowly over the unit cell. The procedure for de- 
riving a theory that is valid close to the Dirac point 
consists in using this representation in the tight-binding 
Hamiltonian and expanding the operators up to a lin- 



J 



ear order in S. In the derivation one uses the fact that 



■£ s e ±iK - s = J2s e = 0- After som e straightfor- 

ward algebra we arrive at 1 Semenofll . Il984f l : 



H ~ -t J dxdyi>\(r) 



3o(l-iV3)/4 
-3a(l-H\/3)/4 

3a(l + i\/3)/4 



-3a(l-i\/3)/4 
*i(r)<7 • V*i(r) + ^\{r)rr* ■ V* 2 (r) ) 

r 



d x - 



3a(-i-\/3)/4 . 

-3a(i-V3)/4 

3a(i-V3)/4 
-3a(-i- V3)/4 



9, 



#i(r) 
* 2 (r) 
(18) 



with Pauli matrices <x = (cr^CTy), <x* = (a x ,— a y ), and 
<3>j = (a|, foj) (i = 1,2). It is clear that the effective 
Hamiltonian lfl8|) is made of two copies of the massless 
Dirac-like Hamiltonian, one holding for p around K and 
other for p around K' . Notice that, in first quantized lan- 
guage, the two-component electron wavefunction, ip(r), 
close to the K point, obeys the 2D Dirac equation: 



(19) 



The wavefunction, in momentum space, for the mo- 
mentum around K has the form: 



^±,K(k) = 



1 



( 



y/2 V ±e 



-i6 k /2 
\9 k /2 



(20) 



for Hk = VFfT • k, where the ± signs correspond to the 
eigenenergies E = ±-upfc, that is, for the 7r and 7r* band, 
respectively, and 9^ is given by (J9j) . The wavefunction 
for the momentum around K' has the form: 



1 / e w */ 2 



(21) 



for Hk' = vpcr* ■ k. Notice that the wavefunctions 
at K and K' are related by time reversal symmetry: 
if we set the origin of coordinates in momentum space 
in the M-point of the BZ (see Figt2j , time reversal be- 
comes equivalent to a reflection along the k x axis, that 
is, (k x ,ky) — > (k x , ~ky). Also note that if the phase 6 is 
rotated by 2ir the wavefunction changes sign indicating 
a phase of ir (in the literature this is commonly called a 
Berry's phase). This change of phase by ir under rotation 
is characteristic of spinors. In fact, the wavefunction is a 
two component spinor. 

A relevant quantity used to characterize the eigenfunc- 
tions is their helicity defined as the projection of the 
momentum operator along the (pseudo)spin direction. 
The quantum mechanical operator for the helicity has 
the form: 



% 1 P 
2 P 



(22) 



It is clear from the definition of h that the states ipK{ r ) 
and ^K'{f) are also eigenstates of h: 



htp K (r) 



±-ip K (r), 



(23) 



and an equivalent equation for 4>K'( r ) with inverted sign. 
Therefore electrons (holes) have a positive (negative) he- 
licity. Equation 11231) implies that <x has its two eigen- 
values either in the direction of or against (JJ.) the 
momentum p. This property says that the states of the 
system close to the Dirac point have well defined chi- 
rality or helicity. Notice that chirality is not defined in 
regards to the real spin of the electron (that has not yet 
appeared in the problem) but to a pseudo-spin variable 
associated with the two components of the wavefunction. 
The helicity values are good quantum numbers as long 
as the Hamiltonian fl8|l is valid. Therefore the existence 
of helicity quantum numbers holds only as an asymptotic 
property, which is well defined close to the Dirac points 
K and K' . Either at larger energies or due to the pres- 
ence of a finite t 1 the helicity stops being a good quantum 
number. 



1. Chira I Tunneling and Klein paradox 

In this section we want to address the scattering of 
chiral electr o ns in t wo dimensions by a squ are barrier 
()Katsnelsonl . l2007bl : iKatsnelson et all l2006h . The one 
dimensional scattering of chiral electrons was discussed 
earlier in the contex t of ca rbon nanotubes l Ando et al\ . 
ll998l : lMcEuen et aZ.I . fl999h 

We start by noticing that by a gauge transformation 
the wavefunction ll20l) can be written as: 



1 



( 



1 



V2 V ±e " 



(24) 



We further assume that the scattering does not mix the 
momenta around K and K' points. In Fig. [6] we depict 
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Figure 6 (Color online) Top: Schematic picture of the scat- 
tering of Dirac electrons by a square potential. Bottom: def- 
inition of the angles <f> and 8 used in the scattering formalism 
in the three regions I, II, and III. 



the scattering process due to the square barrier of width 
D. 

The wavefunction in the different regions can be writ- 
ten in terms of incident and reflected waves. In region I 
we have: 



ipi(r) = — = f \ 

+ 71 ( se ' 



(25) 



with <f> = a,rctaii(ky / k x ) , k x = /c_fcos</>, k y = fcp-sin< 
and kp the Fermi momentum. In region II we have: 



ipn(r) 



a / 1 

; 71 V s ' e% 

+ ±( .} 



v s , e l-0) ) e«-«"+*"> , (26) 
with 9 = arctan(/c y /oa;) and 

& = \/(Vo W/H) ~ (27) 
and finally in region III we have a transmitted wave only: 



4>ui(r) = — 



t I 1 



V2 V se 



oi<l> 



„i(k m x+k y y) 



(28) 



with s = sgn(E) and s' = sgn(E — Vq). The coefficients 
r, a, b and t are determined from the continuity of the 
wavefunction, which implies that the wavefunction has 



to obey the conditions ipi(x = 0, y) = ipn{x = 0,y) and 
ipn(x = D,y) = ipui(x = D,y). Unlike the Schodinger 
equation we only need to match the wavefunction but 
not its derivative. The transmission through the barrier 
is obtained from T(4>) = tt* and has the form: 



cos 2 8 cos 2 i 



cos(Dq x ) cos 4> cos 9} 2 + sin (Dq x ) (1— ss'sin <p sin 0) 2 

(29) 

This expression does not take into account a contribu- 
tion from evanescent waves in region II, which is usually 
negligible, unless the chemical potential in region II is at 
the Dirac energy (see section HV.I|) . 

Notice that T{(f)) = T{—<j>) and for values of Dq x 
satisfying the relation Dq x — rwr, with n an inte- 
ger, the barrier becomes completely transparent since 
T{<p) = 1, independently of the value of <j). Also, for 
normal incidence (</>—>■ and 6 — ► 0) and for any 
value of Dq x one obtains T(0) = 1, and the barrier 
is again totally transparent. This result is a manifes- 
tation of the Klein paradox llCalogeracos and Dombeyl . 
Il999l : lltzykson and Zuberl 12006^ 1 and does not occur for 
non-relativistic electrons. In this latter case and for nor- 
mal incidence, the transmission is always smaller than 
one. In the limit |Vo| > \E\, eq. |29|) has the following 
asymptotic form 



T{4>) 



1 — cos 2 (Dq x ) sin 



(30) 



In Fig. [7] we show the angular dependence of T(0) for 
two different values of the potential Vb; it is clear that 
there are several directions for which the transmission 
is one. Similar calculations w ere done for a graphene 
bilayer I Katsnelson et a?.l . l2006t ) with its most distinctive 
behavior being the absence of tunneling in the forward 
(k y ~ 0) direction. 
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The simplest example of a potential barrier is a square 
potential discussed previously. When intervalley scat- 
tering and the lack of symmetry between sublattices are 
neglected, a potential barrier shows no reflection for elec- 
trons incident in the normal direction ( Katsnelson et all . 
200d ) . Even when the barrier separates regions where the 
Fermi surface is electron like on one side and hole like on 
the other, a normally incident electron continues propa- 
gating as a hole with 100% efficiency. This phenomenon 
is another manifestation of the chirality of the Dirac elec- 
trons within each valley, which prevents backscattering in 
general. The transmission and reflection probabilities of 
electrons at different angles depend on the potential pro- 
file along the barrier. A slowly varying barrier is more 
efficient in reflecting e l ectron s at non-zero incident angles 
(|Cheianov and Fal'kcl 120061 ) . 

Electrons moving through a barrier separating p- and 
n-doped graphene, a p-n junction, are transmitted as 
holes. The relation between the velocity and the mo- 
mentum for a hole is the inverse of that for an elec- 
tron. This implies that, if the momentum parallel to 
the barrier is conserved, the velocity of the quasipar- 
ticle is inverted. When the incident electrons emerge 
from a source, the transmitting holes are focused into 
an image of the source. This behavior is the same as 
that of photon s moving in a medium with negative re- 
flection index i Cheianov et al\ . l2007af ). Similar effects 
can occur in graphene quantum dots, where the in- 
ner and ou ter regions contain electrons and holes, re- 
spectively ( Cserti et all . l2007bh . Note that the fact 
that barriers do not impede the transmission of nor- 
mally incident electrons does not preclude the existence 
of sharp resonances, due to the confinement of elec- 
trons with a finite parallel momentum. This leads to 
the possibility of fabricating quantum dots with po- 
tential barriers (jSilvestrov and Efetovi . l2007h . Finally, 
at half-filling, due to disorder graphe ne can be divided 
in electron and hol e charg e puddles I Katsnelson et all . 
20061 : iMartin et all . 120071) . Transport is determined 
by the trans mission across the p-n j unctions betwee n 
these puddles ijCheianov et aZ.l . l2007bUShklovskiil . l2007i 
There is a rapid progress in the measurement of trans- 
port properties of graphene ribbons with additional top 
ates that play the role of tunable potential barriers 



gates that play the role ot tunable potential barriers 
(lHan et a/.l.l2007l;lHuard et aZll2007tlLemme et aZ.l . l2007t 



lOzvilmaz et all . 120071 : IWilliams et all . l2007h 



A magnetic field and potential fluctuations break both 
inversion symmetry of the lattice and time reversal sym- 
metry. The combination of these effects break also the 
symmetry between the two valleys. The transmission 
coefficient becomes valley dependent, and, in general, 
electrons from different valleys propagate along different 
paths. T his opens the possibility of manipulating the val- 
ley index ( Tworzvdlo et «/.l . l2007t ) (valleytronics) in a way 
similar to the control of the spin in mesoscopic devices 
(spintronics). For large magnetic fields, a p-n junction 
separates regions with different quantized Hall conduc- 
tivities. At the junction, chiral currents can flow at both 



edges ( Abanin and Levitovl . 120071 ) . inducing backscatter- 
ing between the Hall currents at the edges of the sample. 

The scattering of electrons near the Dirac point 
by graphene-superconductor junctions differs from 
Andreev scattering proce ss in normal metals 
( Titov and Beenakked . l2006h . When the distance 
between the Fermi energy and the Dirac energy is 
smaller than the superconducting gap, the supercon- 
ducting interaction hybridizes quasiparticles from one 
band with quasiholes in the other. As in the case of 
scattering at a p-n junction, the trajectories of the 
incoming electron and reflected hole (note that hole here 
is meant as in the BCS theory of superconductivity) are 
different from those in si milar processes in metals with 
only one type of carrier I Bhattachariee and Senguptal . 
120061 : iMaiti and Senguptal . l2007t ). ~ 



2. Confinement and zitterbewegung 

Zitterbewegung, or jittery motion of the wavefunction 
of the Dirac pro blem, occurs when one tries to confine the 
Dirac electrons ( Itzykson and ZuberLl2006t ). Localization 
of a wavepacket leads, due to the Heisenberg principle, to 
uncertainty in the momentum. For a Dirac particle with 
zero rest mass, uncertainty in the momentum translates 
into uncertainty in the energy of the particle as well (this 
should be contrasted with the non-relativistic case where 
the position-momentum uncertainty relation is indepen- 
dent of the energy-time uncertainty relation). Thus, for 
a ultra-relativistic particle, a particle-like state can have 
hole-like states in its time evolution. Consider, for in- 
stance, if one tries to construct a wave packet at some 
time t = 0, and let us assume, for simplicity, that this 
packet has a Gaussian shape of width w with momentum 
close to K: 



V'o(r) 



e -r 2 /(2w 2 ) 



iK-: 



' 7TW 



(31) 



where 4> is spinor composed of positive energy states (as- 
sociated with - 0+,k of ([20|) ). The eigenfunction of the 
Dirac equation can be written in terms of the solution 
(1201) as: 



V(r,t) = 



^ «a,k^a,K(k)t 
J V_J ' a=±l 



(32) 



where a± ; k are Fourier coefficients. We can rewrite f31j) 
in terms of lf32|) by inverse Fourier transform and find 
that: 



a±,k 



V2,/,t 



V4, K (k)0. 



(33) 



Notice that the relative weight of positive energy states 
with respect to negative energy states, |a+/a_|, given 
by (|20l) is one, that is, there are as many positive energy 
states as negative energy states in a wavepacket. Hence, 
these will cause the wavefunction to be delocalized at any 
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time t 7^ 0. Thus, a wave packet of electron-like states 
has hole-like components, a result t hat puzzled many re- 
searc hers in the early days of QED ( Itzvkson and Zubeil 
l200fih . 



E(eV) 




-0.2 - 
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Figure 9 (Color online) (a) Lattice structure of the bilayer 
with the various hopping parameters according to the SWM 
model. The A-sublattices are indi cated by the darker s pheres, 
(b) Brillouin zone. Adapted from iMalard et all . 120071 . 



Figure 8 (Color on line) Energy spectrum (in units of i) for a 
graphene ribbon 600a wide, as a function of the momentum 
k along the ribbon (in units of l/(-\/3a)), in the presence of 
confining potential with Vo = 1 eV, A = 180a. 



Consider the tight-bin ding description (jChen et all . 
l2007al : IPeres et afJ . l2006bh of Sec. III. A I when a potential 
Vi on site Ri is added to the problem: 



H e Villi 



(34) 



graphene layers. The simpl est generalization is a bilayer 
( McCann and Fal'kol . 12006^ . A bilayer is interesting be- 
cause the IQHE shows anomalies , although different from 
those observed in a single layer (jNovoselov et all . l2006t l , 
and also a ga p can open between the c onduction and 
valence band ( McCann and Fal'kol . l2006f ). The bilayer 
structure, with the AB stacking of 3D graphite, is shown 
in FigH 

The tight-binding Hamiltonian for this problem can be 
written as: 



where n, is the local electronic density. For simplicity, we 
assume that the confining potential is ID, that is, that 
Vi vanishes in the bulk but becomes large at the edge 
of the sample. Let us assume a potential that decays 
exponentially away from the edges into the bulk with a 
penetration depth, A. In Fig. [8] we show the electronic 
spectrum for a graphene ribbon of width L = 600a, in 
the presence of a confining potential, 



V{x) = Vo 



e -(x-L/2)/X + e -(L/2-x)/X 



(35) 



where x is the direction of confinement and Vo the 
strength of the potential. One can clearly see that in 
the presence of the confining potential the electron-hole 
symmetry is broken and, for Vq > 0, the hole part of 
the spectrum is strongly distorted. In particular, for k 
close to the Dirac point, we see that the hole disper- 
sion is given by: i?„ jCr= _i(fc) sw — 7„fc 2 — ( n k 4 where n 
is a positive integer, and 7„ < (7„ > 0) for n < N* 
(n > N*). Hence, at n = N* the hole effective mass di- 
verges (7at» = 0) and, by tuning the chemical potential, 
jji, via a back gate, to the hole region of the spectrum 
(/i < 0) one should be able to observe an anomaly in the 
Shubnikov-de Haas (SdH) oscillations. This is how zit- 
terbewegung could manifest itself in magnetotransport. 



C. Bilayer graphene: tight-binding approach 



n 



t.b. 



-70 ( a L,i,a b m,j,a + h- c -) 

<i,3> 
rn.a 



73 



74X^( 6 tj> &2 »^> + h - c -)' 



(36) 



where a m i jCr (b m ^ a ) annihilates an electron with spin 
a, on sublattice A (B), in plane m = 1,2, at site R^. 
Here we use the graphite nomenclature for the hop- 
ping parameters: 70 = t is the in-plane h opping energy 
and 7 1 (71 = tj_ sa 0.4 eV in graphite ijBrandt et all . 
Il988l : iDresselhaus and Dresselhausl . l2Q02h ) is the hop- 
ping energy between atom Ai and atom A2 (see Fig. [9t , 
and 73 (73 ~ 0.3 eV in graphite ( Brandt et all . Il988t 
IDresselhaus and Dresselhausl . l2002t )) is the hopping en- 
ergy between atom Ai (A 2 ) and atom B 2 (Bi ), and 
74 (74 ^ —0-04 eV in gra p hite J Brandt et all . 



IDresselhaus and Dresselhausl . 120021 )) that connects Bi 
and B 2 . 

In the continuum limit, by expanding the momentum 
close to the K point in the BZ, the Hamiltonian reads, 



The tight-binding model developed for graphite can 
be easily extended to stacks with a finite number of 



k 



(37) 
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Figure 10 (Color online) Band structure for bilayer graphene 
for V = and 73 = 0. 



Figure 11 (Color online) Band structure for bilayer graphene 
for V and 73 = 0. 



where (ignoring 74 for the time being) : 

/ — V vpk 3^ak* 

( v F k* -V 7l 

A_ 71 V v F k 
\373afc v F k* V 



(38) 



where k = k x + ik y is a complex number, and we have 
added V which is here half the shift in electro-chemical 
potential between the two layers (this term will appear 
if a potential bias is applied between the layers) , and 



4 = M(k) >£ 4(k),6t(k),4(k) 



(39) 



is a four component spinor. 

If V = and ^,vpk <C 71, one can eliminate the 
high energy states perturbatively and write an effective 
Hamiltonian: 



Hk = 







71 



*4r 
71 



373 a A; 



I- 37 3 a/c* 




(40) 



The hopping 74 leads to a fc dependent coupling between 
the sublattices or a small renormalization of 71. The 
same role is played by the inequivalence between sublat- 
tices within a layer. 

For 73 = 0, (|40|) gives two parabolic bands, ek,± ~ 
/ijL which touch at e = (as shown in Fig fT0|) . 
The spectrum is electron-hole symmetric. There are two 
additional bands which start at ±ij_. Within this approx- 
imation, the bilayer is metallic, with a constant density of 
states. The term 73 changes qualitatively the spectrum 
at low energies since it introduces a trigonal distortion, 
or warping, of the bands (notice that this trigonal dis- 
tortion, unlike the one introduced by large momentum 
in ((8j), occurs at low energies). The electron-hole sym- 
metry is preserved but, instead of two bands touching at 
k = 0, we obtain three sets of Dirac-like linear bands. 
One Dirac point is at e = and k = 0, while the three 
other Dirac points, also at e = 0, lie at three equivalent 
points with a finite momentum. The stability of points 
where bands touch can be understood using topological 



arguments ( Manes et aZ.l . l2007h . The winding number of 
a closed curve in the plane around a given point is an 
integer representing the total number of times that the 
curve travels counterclockwise around the point so that 
the wavefunction remains unaltered. The winding num- 
ber of the point where the two parabolic bands come 
together for 73 = has winding number +2. The trigo- 
nal warping term, 73, splits it into a Dirac point at k = 
and winding number —1, and three Dirac points at k 7^ 
and winding numbers +1. An in-plane magnetic field, or 
a small rotation of one layer with respect to the other 
splits the 73 = degeneracy into two Dirac points with 
winding number +1. 

The term V in f38|) breaks the equivalence of the two 
layers, or, alternatively, inversion symmetry. In this case, 
the dispersion relation becomes: 



V 2 



2 7 2 

Vpk 



ti/2 



± 



sJw 2 v 2 F k 2 



t 2 v 2 F k 2 



(41) 



given rise to the dispersion shown in Fig. [TTJ and to the 
opening of a gap close, but not directly at, the K point. 
For small momenta, and V <C t, the energy of the con- 
duction band can be expanded: 



e k ~ V 



(2Vv 2 F k 2 )/t A 



+ (v F k i )/(2t 2 ± V). 



(42) 



The dispersion for the valence band can be obtained by 
replacing by — efe. The bilayer has a gap at k 2 sa 
(2V 2 )/v F . Notice, therefore, that the gap in the biased 
bilayer system depends on the ap plied bias an d hence 
can be measured experimentally d Castro et all l2007at 
iMcCannl . 120061 : iMcCann and Fal'kol . l2006t) . The ability 
to open a gap makes bilayer graphene most interesting 
for technological applications. 



D. Epitaxial graphene 

It has been known for a long time that monolay- 
ers of graphene could be grown epitaxially on metal 
surfaces by using catal ytic decomposition of hvdrocar- 
bons or carbon oxide 1 Campagnoli and Tosattil . Il989t 



12 



Eizcnberg and BlakelvL Il97£ 



199/ 
200f 



Shelton et all 



Oshima and Nagashima, 
ISinitsyna and Yaminsky . 



When such surfaces are heated, oxygen or hy- 
drogen desorbs, and the carbon atoms form a graphene 
monolayer. The resulting graphene structures could 
reach sizes up to a micrometer, with few defects and were 
characterized by differ ent surface-science te chniques and 
local scanning probes IjHimpsel et all 1 1982h . For exam- 
ple, graphene grown on ruthenium has zigzag edges and 
also ripples associated with a (10 x 10) reconstruction 
([Vazquez de Parga et all l2007l l. 

Graphene can also be formed on the surface of 
SiC. Upon heating, the silicon from the top layers 



the surface dBerger et al. 2004: Bommel et al. 


. 1975 


Coev et al. 20021: ] 


7 orbeaux et all Il998|: Hass et al. 


2007a: Ide Heer et al. 


. 20071: Rollings et all l2005h. The 



number of layers can be controlled by limiting time or 
temperature of the heating treatment. The quality and 
the number of layers in th e samples depends o n the SiC 
face used for their growth I de Heer et all [20071 ) (the car- 
bon terminated surface produces few layers but with 
a low mobility, whereas the silicon terminated surface 
produces several layers but with higher mobility). Epi- 
taxially grown multilayers e xhibit SdH o s cillati ons with 
a Berry phase shift of n ( Berger et all l2006h . which 
is the same as the phase shift for Dirac fermions ob- 
served in a single layer as well as for some subbands 
present in multilayer graphene (see fu rther) and graphite 
( Luk'yanchuk and Kopelevichl . 12004^ . The carbon layer 
directly on top of the substrate is expecte d to be strongly 
bond ed to it, and it shows no it bands ijVarchon et all 
12007( 1. The next layer shows a (6^3 x 6^3) reconstruc- 
tion due to the substrate, and has graphene properties. 
An alternate route to produce few laye rs graphene is 
based on synthesis from nanodiamonds ( Affoune et all 
12001( 1 . 

Angle resolved photo-emission experiments (ARPES) 
show that epitaxial graphene grown on SiC has 
l inearly disper s ing quasiparticles (Dirac fermions) 
dBostwick et all l2007rl lOhta et all l2007t IZhou et all 
l2006bh . in agreement with the theoretical expectation. 
Nevertheless, these experiments show that the electronic 
properties can change locally in space indicating a cer- 
tain degree of inho mogeneity due to the growth method 
((Zhou et aZ.I . l2007f ). Similar inhomogeneities due to dis- 
order in the c-axis o rientation o f graph ene planes is 
observed in graphite I Zhou et~al\ . l2006a( l. Moreover, 
graphene grown this way is heavily doped due to the 
charge transfer from the substrate to the graphene layer 
(with the chemical potential well above the Dirac point) 
and therefore all samples have stro ng metallic charac- 
ter with large electr onic mobilities ijBerger et all 120061 : 
Ide Heer et all l2007t l . There is also evidence for strong 
interaction between a substrate and the graphene layer 
l eading to the ap pearance of gaps at the Dirac point 
( Zhou et~ai\ . l2007( l. Indeed, gaps can be generated by 
the breaking of the sublattice symmetry and, as in the 




Figure 12 (Color online) Sketch of the three inequivalent ori- 
entations of graphene layers with respect to each other. 



case of other carbon based systems such as polyacethy- 
lene (ISu et all Il979l. Il980h. it can lead to soliton-like 



excitations ( Hou et alT 20071 : Ijackiw and Rebbl Il976f l . 



Multilayer gra phene grown on SiC have also been s t udied 
with ARPES jBostwick eTal\ . l2007at : lOhta et all 120071 
I200(t1 1 and the results s eem to agree quite well with ban d 
structure calculations ( Mattausch and Pankratovl . l2007( l. 
Spectroscopy measurements also show the transitio ns as- 
sociated with Landau levels ( Sadowski et all feoOtjf ) . and 
weak localization effects at low magnetic fields, also ex- 
pected for Dirac fermions (|Wu et all 20071). Local probes 
reveal a rich structure of terraces dMallet et all l2007h 



and interference patterns due to def ects at or below the 
graphene layers ( Rutter et all l2007( l . 



E. Graphene stacks 

In stacks with more than one graphene layer, two con- 
secutive layers are normally oriented in such a way that 
the atoms in one of the two sublattices, A n , of the honey- 
comb structure of one layer are directly above one half of 
the atoms in the neighboring layer, sublattice A n ±i. The 
second set of atoms in one layer sits on top of the (empty) 
center of an hexagon in the other layer. The short- 
est distance between carbon atoms in different layers is 
dA n A n±1 = c = 3.4A. The next distance is c?a„s„ ± i = 
Vc 2 + a 2 . This is the most common arrangement of near- 
est neighbor layers observed in nature, although a stack- 
ing order in which all atoms in one layer occupy posi- 
tions directly above the atoms in the neighboring layers 
(hexagonal sta c king) has been considered theoretically 
( Charlier et aD, Il99lh and appears in graphite interca- 
lated compounds ([Dresselhaus and Dresselhausl . 120021 ) . 

The relative position of two neighboring layers allows 
for two different orientations of the third layer. If we la- 
bel the positions of the two first atoms as 1 and 2, the 
third layer can be of type 1, leading to the sequence 121, 
or it can fill a third position different from 1 and 2 (see 
Fig. [12]), labeled 3. There are no more inequivalent po- 
sitions where a new layer can be placed, so that thicker 
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Figure 13 (Color online) Electronic bands of graphene multi- 
layers: top left: biased bilayer; top right: trilayer with Bernal 
stacking; bottom left: trilayer with orthorhombic stacking; 
bottom right: stack with four layers where the top and bot- 
tom layers are shifted in energy with respect to the two middle 
layers by +0.1 eV. 



stacks can be described in terms of these three orien- 
tations. In the most common version of bulk graphite 
the stacking order is 1212 •• • (Bernal stacking). Re- 
gions with the stacking 123123 • • • (rhombohedral stack- 
i ng) have also been observed in d ifferent types of graphite 
( Baconl . [l950l : iGasparouxL Il967h . Finally, samples with 
no discernible stacking order (turbostratic graphite) are 
also commonly reported. 

Beyond two layers, the stack ordering can be arbi- 
trarily complex. Simple analytical expressions for the 
electronic bands can be obtained for perfect Bernal ( 
1212 ■ ■ ■ ) and rho mbohedral ( 123123 • • • ) stacking 
( Guinea et iCven if we consider one interlayer 

hopping, t± = 71, the two stacking orders show rather 
different band structures near e = 0. A Bernal stack 
with N layers, TV even, has N/2 electron like and N/2 
hole like parabolic subbands touching at e = 0. When 
N is odd, an additional subband with linear (Dirac) dis- 
persion emerges. Rhombohedral systems have only two 
subbands that touch at e = 0. These subbands disperse 
as k N , and become surface states localized at the top 
and bottom layer when TV — > oo. In this limit, the re- 
maining 2N — 2 subbands of a rhombohedral stack be- 
come Dirac like, with the same Fermi velocity as a sin- 
gle graphene layer. The subband structure of a tri-layer 
with the Bernal stacking includes two touching parabolic 
bands, and one with Dirac dispersion, combining the fea- 
tures of bilayer and monolayer graphene. 

The low energy bands have different weights on the two 
sublattices of each graphene layer. The states at a site 
directly coupled to the neighboring planes are pushed to 
energies e ~ ±t±. The bands near e = are localized 
mostly at the sites without neighbors in the next layers. 
For the Bernal stacking, this feature implies that the den- 
sity of states at e = at sites without nearest neighbors in 
the contiguous layers is finite, while it vanishes linearly 
at the other sites. In stacks with rhombohedral stack- 



ing, all sites have one neighbor in anot her plane, and 
the d ensity of states vanishes at e = ( Guinea et all . 
20061 ). This result is consistent with the well known fact 
that only one of the two sublattices at a graphite surface 
can be resolved by scan ning tunneling microscopy (STM) 



(jTomanek et all . Il987l ) 



As in the case of a bilayer, an inhomogeneous charge 
distribution can change the electrostatic potential in the 
different layers. For more than two layers, this break- 
ing of the equivalence between layers can take place even 
in the absence of an applied electric field. It is interest- 
ing to note that a gap can open in a stack with Bernal 
ordering and four layers, if the electronic charge at the 
two surface layers is different from that at the two in- 
ner ones. Systems with a higher number of layers do not 
show a gap, even in the presence of charge inhomogene- 
ity. Four representative examples are shown in Fig. \TE[ 
The band structure analyzed here will be modified by 
the inclusion of the trigonal warping term, 73. Experi- 
mental studies of graphene stacks have showed that, with 
increasing number of layers, the system becomes increas- 
ingly metallic (concentration of charge carriers at zero en- 
ergy gradually increases), and the re appear severa l types 
of electron-and-hole-li ke carries ( Morozov et all . 120051 : 
iNovoselov et al I l2004l h An inhomogeneous charge dis- 
tribution between layers becomes very important in this 
case, leading to 2D electron and hole systems that occupy 
only a few graphene layers near the surface and can com- 
pletely dominate tran sport properties of graphene stacks 
(jMorozov eFaIl . l2005h . 

The degeneracies of the ba nds at e = can be studied 
using topological arguments (jMafies et all . l2007t ) . Multi- 
layers with an even number of layers and Bernal stacking 
have inversion symmetry, leading to degeneracies with 
winding number +2, as in the case of a bilayer. The 
trigonal lattice symmetry implies that these points can 
lead, at most, to four Dirac points. In stacks with an odd 
number of layers, these degeneracies can be completely 
removed. The winding number of the degeneracies found 
in stacks with TV layers and orthorhombic ordering is ±N. 
The inclusion of trigonal warping terms will lead to the 
existence of many weaker degeneracies near e = 0. 

Furthermore, it is well known that in graphite the 
planes can be rotated relative each other giving rise to 
Moire patterns that are observ ed in STM of graphite sur- 
faces ( Rong and Kuiperl . Il993l ) . The graphene layers can 
be rotated relative to each other due to the weak coupling 
between planes that allows for the presence of many dif- 
ferent orientational states that are quasidegenerate in en- 
ergy. For certain angles the graphene layers become com- 
mensurate with each other leading to a lowering of the 
electronic energy. Such phenomenon is quite similar to 
the commensurate-incommensurate transitions observed 
in certain charge density wav e systems or adsorption of 
gases on graphite ( Bakl . Tl982h . This kind of dependence 
of the electronic structure on the relative rotation angle 
between graphene lay ers leads to what is cal led super- 
lubricity in graphite (jPienwiebel et all . \2004) . namely, 
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the vanishing of the friction between layers as a func- 
tion of the angle of rotation. In the case of bilayer 
graphene, a rotation by a small commensurate angle leads 
to the effective decoupling between layers and the recov- 
ery of the linear Dirac spectrum of the single layer albeit 
with a modification on the va lue of the Fermi velocity 
( Lopes dos Santos et q/.I . [20Q7h . 



1. Electronic structure of bulk graphite 

The tight-binding description of graphene described 
earlier can be extended to systems with an infinite num- 
ber of layers. The coupling between layers leads to hop- 
ping terms between tt orbitals in different layers, lead- 
i ng to the so called Slqnczew ski-Weiss-McClure model 
( Slonczewski and Weissl . fl958l ). This model describes the 
band structure of bulk graphite with the Bernal stacking 
order in terms of seven parameters, 70,71,72,73,74,75 
and A. The parameter 70 describes the hopping within 
each layer, and it has been considered previously. The 
coupling between orbitals in atoms that are nearest 
neighbors in successive layers is 71, which we called t± 
earlier. The parameters 73 and 74 describe the hop- 
ping between orbitals at next nearest neighbors in suc- 
cessive layers and were discussed in the case of the bi- 
layer. The coupling between orbitals at next nearest 
neighbor layers are 72 and 75. Finally, A is an on site 
energy which reflects the inequivalence between the two 
sublattices in each graphene layer once the presence of 
neighboring layers is taken into account. The values of 
these parameters, and their dependence with pressure, 
or, equivalently, t he interatomic distances, have been ex- 
tensively studied llBrandt afJ.ll988HDillon _et a*J . |l977 : 



Dresselhaus and Mavroidesl . I1964J : iMcCfurel . 1 19571 11964 : 
Nozieresl . Il958t ISoule et all ll964T ) . A representative set 



al. 

37, 



of values is shown in Table [T]. It is unknown, however, 
how these parameters may vary in graphene stacks with 
a small number of layers. 

The unit cell of graphite with Bernal stacking includes 
two layers, and two atoms within each layer. The tight- 
binding Hamiltonian described previously can be repre- 
sented as a 4 x 4 matrix. In the continuum limit, the 
two inequivalent corners of the BZ can be treated sepa- 
rately, and the in plane terms can be described by the 
Dirac equation. The next terms in importance for the 
low energy electronic spectrum are the nearest neighbor 
couplings 71 and 73 . The influence of the parameter 74 
on the low energy bands is much smaller, as discussed 
below. Finally, the fine details of the spectrum of bulk 
graphite are determined by A, which breaks the electron- 
hole symmetry of the bands preserved by 70,71 and 73. 
It is usually assumed to be much smaller than the other 
terms. 

We label the two atoms from the unit cell in one layer 
as 1 and 2, and 3 and 4 correspond to the second layer. 
Atoms 2 and 3 are directly on top of each other. Then, 
the matrix elements of the Hamiltonian can be written 



70 


3.16 eV 


71 


0.39 eV 


72 


-0.020 eV 


73 


0.315 eV 


74 


-0.044 eV 


75 


0.038 eV 


A 


-0.008 eV 



Table I Band structur e pa rameters of graphite 
ijDresselhaus and "Dr esselhausl 120021 ') . 



as: 



H 12 



H 22 



H23 = 



H 24 
H 34 



272 cos(27rfc z /c) 
v F (k x + iky) 

(1 + e ik * c ) (k x + iky) 

(1 + C lk ' C ) (k X ~ iky) 

A + 27 5 cos(27rfc z /c) 
7i(l + e^ c ) 

(1 + e lk > c ) (k x + iky) 



374a 
2 

373 a 



374a 



A + 275 cos(27rfc z /c) 
v F (k x + iky) 
272 cos(2irk z / c) 



(43) 



where c is the lattice constant in the out of plane direc- 
tion, equal to twice the interlayer spacing. The matrix 
elements of H can be obtained by replacing k x by — k x 
(other conventions for the unit cell and the orientation 
of the lattice lead to differen t phase s ) . Recent ARPES 
experiment s llBostwick et all l2007bl ; Ohta et Ml, [200i 
IZhou et all l2006allch performed in epitaxially grown 
graphene stacks ( Berger et all . l2004h confirm the main 



features of this model, formulated m ainly on the ba- 
sis of Fermi su rface measurements ( McClurel . Il957t 
ISoule et a/.l . ll964| l. The electronic spectrum of the model 
can also be c alcula ted in a magnetic field ( de Gennell . 
Il964l : iNakad I 1976h . and the res ults are also con s istent 
with STM on g raphi t e surfaces llKobavashi et al. L 200!: ; 



Li and Andrei . l2007l : iMatsui et all 120051 : iNiimi et al 



2006), epitaxially grown graphene stacks (|Mallet et al 
20071 ). and with opt ical measurements in the infrared 
range (|Li et a/J . [2006h . 



F. Surface states in graphene 

So far, we have discussed the basic bulk proper- 
ties of graphene. Nevertheless, graphene has very in- 
teresting surface (edge) states that do not occur in 
other systems. A semi- infinite graphene sheet with a 
zigzag edge has a band of z ero energy states lo c alized 
at the surface i Fuiita et all Il996t iNakada et all Il996t 
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Figure 14 (Color online) Ribbon geometry with zigzag edges. 



IWakabavashi et all . Il999f l . In section III. HI we will dis- 
cuss the existence of edge states using the Dirac equa- 
tion. Here will discuss the same problem using the tight- 
binding Hamiltonian. To see why these edge states exist 
we consider the ribbon geometry with zigzag edges shown 
in Fig. HH The ribbon width is such that it has N unit 
cells in the transverse cross section (y direction). We will 
assume that the ribbon has infinite length in the longi- 
tudinal direction (x direction). 

Let us rewrite ([5]), with t' — 0, in terms of the integer 
indices m and n, introduced in Fig. HH and labeling the 
unit cells: 

H = —t [a£.(m, n)& ff (m, n) + a\(m, n)b a (m — 1, n) 

m,n,er 

+al(m, n)6 CT (m, n - 1) +h.c.]. (44) 

Given that the ribbon is infinite in the ai direction one 
can introduce a Fourier decomposition of the operators 
leading to 

H = -t [ ^$}4(ft,n)M*>») + e ika al(k,n)b a (k,n) 



+a'l(k, n)ba(k, n — 1) + h.c. 



where cj.(fc, n) |0) = |c, a, k, n), and c = a,b. The one- 
particle Hamiltonian can be written as: 

H lp = -t J dkY^i^ + e ika )\a,k,72,a) (b,k,n,a\ 

+ \a, k, n, a) (b, k,n-l,a\ + h.c.]. (46) 

The solution of the Schrodinger equation, H lp |/i, fc, a) = 
Efi,k | A*, A;, tx) , can be generally expressed as: 

fc, a) = ^^[a(fc, n) |a, k, n, a) + /3(fc, n) \b, k, n, tr)], 

n 

(47) 

where the coefficients a and f3 satisfy the following equa- 
tions: 

E^ k a(k,n) = -t[(l + e lka )(3(k,n)+f3(k,n~l)},(48) 
E^ k /3(k,n) = -t[(l + e~ lka )a(k,n)+a{k,n+lj^9) 



As the ribbon has a finite width we have to be careful 
with the boundary conditions. Since the ribbon only ex- 
ists between n = and n = N — 1 at the boundary 
Eqs. (@HJ| and (02} read: 



E^ k a{k,Q) 
E^ k p(k,N-l) 



-t(l + e ifca )/3(fc,0), (50) 
-t{l + e-' lka )a(k,N - 1). (51) 



The surface (edge) states are solutions of Eqs. (|48H5ip 
with E^k = 0: 







= (1 



(1 + e lka )p(k,n) +P(k 7 n- 1) 



l )a(k, n) + a(k, n + 1) . 



= /?(fc,0)> 

= a(k,N-i). 

Equations |52|) and (|55|l are easily solved giving: 



«(*,») = [-2cos(/fca/2)]"e^ n a(fc,0), 



(52) 
(53) 
(54) 
(55) 



(56) 



/3(fc,n) 



-2cos(fca/2)] 



JV-l-n, 



'^-^^(fc,^-^) 



Let us consider, for simplicity, a semi-infinite system 
with a single edge. We must require the convergence con- 
dition |— 2 cos(fea/2)| < 1, in (|57|) because otherwise the 
wavefunction would diverge in the semi-infinite graphene 
sheet. Therefore, the semi-infinite system has edge states 
for ka in the region 2ir/3 < ka < 47r/3, which corre- 
sponds to 1/3 of the possible momenta. Note that the 
amplitudes of the edge states are given by, 



\a(k,n)\ 
\(3(k,n)\ 



X(k) 



X(k) 



,-n/X(k) 



-(JV-l-n)/A(fc) 



(45) where the penetration length is given by: 



X(k) = -l/ln|2cos(fca/2)|. 



(58) 



(59) 



(60) 



It is easily seen that the penetration length diverges when 
ka approaches the limits of the region ]2n/3, 47r/3[. 

Although the boundary conditions defined by Eqs. lf54j) 
and ([55| are satisfied for solutions ([56| and |57|) in the 
semi-infinite system, they are not in the ribbon geome- 
try. In fact, Eqs. f58|) and f59|) are eigenstates only in 
the semi-infinite system. In the graphene ribbon the two 
edge states, which come from both sides of the edge, will 
overlap with each other. The bonding and anti-bonding 
states formed by the two edge states will th en be the rib- 
bon eigenstates I Wakabavashi et all . ll999h (note that at 
zero energy there are no other states with which the edge 
states could hybridize). As bonding and anti-bonding 
states result in a gap in energy the zero energy flat bands 
of edge states will become slightly dispersive, depending 
on the ribbon width N. The overlap between the two 
edge states is larger as ka approaches 27r/3 and 47r/3. 
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Figure 15 (Color online) Sketch of a zigzag terminati on of a 
graphene bilayer. As discussed in ijCastro et a/.Ll2007bh . there 
is a band of surface states completely localized in the bottom 
layer, and another surface band which alternates between the 
two. 



This means that deviations from zero energy flatness will 
be stronger near these points. 

Edge states in graphene nanoribbons, just as the case 
of carbon nanotubes, are predicted to be Luttinger liq- 
uids, that is, interacting one-d imensional electron sys- 
tems I Castro Neto et a/Xl2006bh . Hence, clean nanorib- 
bons must hav e ID square root sing ularities in their den- 
sity of states ( Nakada et all . Il996h that can be probed 



by Raman spectroscopy. Disorder may smooth out these 
singularities, however. In the presence of a magnetic 
field, when the bulk states are gapped, the edge states 
are responsibl e for t he transport of spin and charge 
(lAbanin et all I2006L l2007al : lAbanin and Levitovl . 120071 : 
lAbanin et aZTl2007bh . 



G. Surface states in graphene stacks 

Single layer graphene can be considered a zero gap 
semiconductor, which leads to the extensively studied 
possibility of gap states, at e = 0, as discussed in the 
previous section. The most studied s uch states are those 
localized near a graphene zigzag edge i Fuiita et «/J . ll996l : 
IWakavabashi and Sigristl. l2000h . It can be shown analyt- 
ically (jCastro et all . l2007bT l that a bilayer zigzag edge, 
like that shown in Fig. HH analyzed within the nearest 
neighbor tight-binding approximation described before, 
has two bands of localized states, one completely local- 
ized in the top layer and indistinguishable from similar 
states in single layer graphene, and another band which 
alternates between the two layers. These states, as they 
lie at e = 0, have finite amplitudes on one half of the sites 
only. 

These bands, as in single layer graphene, occupy one 
third of the BZ of a stripe bounded by zigzag edges. They 
become dispersive in a biased bilayer. As graphite can 
be described in terms of effective bilayer systems, one 
for each value of the perpendicular momentum, fc z , bulk 
graphite with a zigzag termination should show one sur- 
face band per layer. 




zigzag edge 



Figure 16 (Color online) A piece of a honeycomb lattice dis- 
playing both zigzag and armchair edges. 



H. The spectrum of graphene nanoribbons 

The spectrum of graphene nanoribbons depend very 
much on the nature of t he ir edges - zigzag or armchair 
(|Brev and Fertigl . l2006al |bl: iNakada et a/1 . Il996h . In Fig. 
[TBI we show a honeycomb lattice having zigzag edges 
along the x direction and armchair edges along the y 
direction. If we choose the ribbon to be infinite in the x 
direction we produce a graphene nanoribbon with zigzag 
edges; conversely choosing the ribbon to be macroscop- 
ically large along the y but finite in the x direction we 
produce a graphene nanoribbon with armchair edges. 

In Fig. [17] we show the fourteen energy levels, calcu- 
lated in the tight-binding approximation, closest to zero 
energy for a nanoribbon with zigzag and armchair edges 
and of width A*" = 200 unit cells. We can see that they 
are both metallic, and that the zigzag ribbon presents a 
band of zero energy modes that is absent in the armchair 
case. This band at zero energy is the surface states liv- 
ing near the edge of the graphene ribbon. More detailed 
oh initio calculations of the spectra of graphene nanorib- 
bons show that interac tion effects can lead to electronic 
gaps l|Son et all . [2006bh and magnetic states close to the 
graphene edges, independ ent of their nature i Son et all . 
l2006at lYang et all . l2007al |bl) . 

From the experimental point of view, however, 
graphene nanoribbons currently have a high degree 
of roughness at the edges. Such edge disorder 
can change significantly the properties of edge states 
( Areshkin and Whitd . l2007l : lGunivcke et artl2007t h lead- 
ing to Anderson localiza tion, and anom al ies in 
the quantum Hall effect ( Castro Neto et all . l2006bt 
iMartin and Blanterl . l200?t l as well as Coulomb block- 
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and around the K' as: 



1 2 3 4 5 
momentum ka 



2.5 3 3.5 
momentum ka 



Figure 17 (Color online) Left: Energy spectrum, as calcu- 
lated from the tight-binding equations, for a nanoribbon with 
armchair(top) and zigzag(bottom) edges. The width of the 
nanoribbon is N = 200 unit cells. Only fourteen eigenstates 
are depicted. Right: Zoom of the low energy states shown on 
the right. 



ade effects 1 Sols et all . l2007f l. Such effects have already 
been observe d in lithog r aphica l ly engineered g r aphen e 
nanoribbons (jHan et all 120071 : lOzvilmaz et all 120071 ). 



Furthermore, the problem of edge passivation by hydro- 
gen or other elements is not clearly understood experi- 
mentally at this time. Passivation can be modeled in the 
tight-bin ding approach b y modifications of the hopping 
energies (|Novikovll27i(T73f or via addit i onal p hases in the 
boundary conditions ( Kane and Meld . Il997h . Theoreti- 
cal modeling of edge passivation indicate that those have 
a strong effect on the el ectronic properties at the edge 
of gra phene nanoribbons i Barone et all . [20061 : iHod et all 
l2007h . 

In what follows we derive the spectrum for both 
zigzag and armchair edges directly from the Dirac equa- 
tion. This was originally done both wi th and without a 
magn etic field ( Brey and Fertigl . l2006al lbl: iNakada et all . 
1996). 



1. Zigzag nanoribbons 

In the geometry of Fig.[l6]the unit cell vectors are a\ = 
ao(l, 0) and a-i = (l/2, s/3/2) , which generate the unit 
vectors of the BZ given by 6i = 47r/(a \/3) (y/3/2, -1/2) 
and £>2 = 47r/(aov / 3)(0, 1). From these two vectors 
we find two inequivalent Dirac points given by K = 
(47r/3a ,0) = (K, 0) and K' = (-47r/3a ,0) = (-K,0), 
with ao = \/3a. The Dirac Hamiltonian around the Dirac 
point K reads in momentum space: 



H 



K 



Vp 



p x - iPy 

Px + iPy 



(61) 



H K , = v F 



P X + IPy 

Px ~ iPy 



(62) 



The wavefunction, in real space, for the sublattice A is 
given by: 



* A (r) = e iK - r Mr) + e iJC '" r ^(r) , 
and for sublattice B is given by 



*u(r) - e^ r Mr) + e lK '- r ^' B {r) , (64) 

where ipA and ipB are the components of the spinor wave- 
function of Hamiltonian lj6lj) and ip' A and ip' B have iden- 
tical meaning but relatively to ([62|l. Let us assume that 
the edges of the nanoribbons are parallel to the a;— axis. 
In this case, the translational symmetry guarantees that 
the spinor wavefunction can be written as: 



ip{r) = e 



(65) 



and a similar equation for the spinor of Hamiltonian (|62|l . 
For zigzag edges the boundary conditions at the edge of 
the ribbon (located at y — and y = L, where L is the 
ribbon width) are: 



* A (y = L)=0, 
leading to: 



*B(y = 0) = 0, 



(66) 



= e iKx e ik * x <p A (L) + e- iKx e ik * x <j> A (L), (67) 

q „iKx „ik x x x _ ( f\\ i „ — iKx „ik x x jj 



AKx jik^x , 



e^^' B (0). (68) 



The boundary conditions |67|) and ([68)l are satisfied for 
any x by the choice: 

^(i)=^(L) = ^ s (0) = ^(0)=0. (69) 

We need now to find out the form of the envelope func- 
tions. The eigenfunction around the point K has the 
form: 



o K-d y \( My) \ = J My) 
kx + d y o j\ My) J V Mv) 



(70) 



with e = e/vp and e the energy eigenvalue. The eigen- 
problem can be written as two linear differential equa- 
tions of the form: 



(71) 



(k x - dy)4> B = e<j>A ■, 
(k x + dy)<p A = e(j> B - 



Applying the operator [k x + d y ) to the first of Eqs. (71 
leads to: 



{-d 2 y + kl)^ B = i^ B: 



(72) 
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with 4>a given by: 



4>a = \{k x - d y )4>B 



The solution of i|72]) has the form: 

(j) B = Ae zy + Be~ zy , 



(73) 



(74) 



leading to an eigenenergy e 2 — k% — z 2 . The boundary 
conditions for a zigzag edge require that <j) A {y = L) =0 
and 4*b(jj = 0) =0, leading to: 

f 4>b{v = 0) = 0^ A + B = Q, 

[ <f> A {y = L) = (fc K — z)Ae zL + (fc x + z)Ber zL = 

(75) 

which leads to an eigenvalue equation of the form: 



-2zL 



k x + z 



(76) 



Equation (|76j) has real solutions for z, whenever k x is pos- 
itive; these solutions correspond to surface waves (edge 
states) existing near the edge of the graphene ribbon. In 
section III.FI we discussed these states from the point of 
view of the tight-binding model. In addition to real solu- 
tions for z, l(76"j) also supports complex ones, of the form 
z = ik n , leading to: 



kn 



tan(fc n L) 



(77) 



The solutions of f77j) correspond to confined modes in 
the graphene ribbon. 

If we apply the same procedure to the Dirac equation 
around the Dirac point K ' we obtain a different eigen- 
value equation given by: 



-2zL 



k x + z 



(78) 



This equation supports real solutions for z if k x is neg- 
ative. Therefore we have edge states for negative values 
k x , with momentum around K 1 . As in the case of K, 
the system also supports confined modes, given by: 



k x — 



k n 



tan(fc„L) 



(79) 



One should note that the eigenvalue equations for K' are 
obtained from those for K by inversion, k x — + — k x - 

We finally notice that the edge states for zigzag 
nanoribbons are dispersionless (localized in real space) 
when t' = 0. When electron-hole symmetry is broken 
(f 7^ 0) these s tates become dispersive w ith a Fermi ve- 
locity v e w t'a (|Castro Neto et all l2006bh . 



2. Armchair nanoribbons 

Let us now consider an armchair nanoribbon with arm- 
chair edges along the y direction. The boundary condi- 
tions at the edges of the ribbon (located at x = and 



x = L, where L is the width of the ribbon): 

V A (x = 0) = V B {x = 0) = y A (x = L) = V B (x = L) = . 

(80) 

Translational symmetry guarantees that the spinor wave- 
function of Hamiltonian il6ll) can be written as: 



(81) 



\ <I>b(x) 



and a similar equation for the spinor of the Hamiltonian 
<[62|l . The boundary conditions have the form: 



= e lk y y (j) A (0) + e lk y y (t>' A {0) 



(82) 

= e^^BfOj + ^AfO), (83) 
= e iKL e ik y y <t> A (L) + e- iKL e ik v v (f> A (L), (84) 
= e lKL e lk y y (t)B{L) + e- lKL e tk y y (t)' B {L), (85) 



and are satisfied for any y if: 

^(0) + ^(0) = ; 

and 



iKL 



ML) 



-iKL 1 1 



(86) 



(87) 



with [i = A, B. It is clear that these boundary conditions 
mix states from the two Dirac points. Now we must find 
the form of the envelope functions obeying the boundary 
conditions (|86j) and (|87|) . As before, the functions cj>B 
and 4>' B obey the second order differential equation f72|) 
(with y replaced by x) and the function <p A and (j>' A are 
determined from <[73|l . The solutions of ([72| have the 
form: 



: >B 



Ae ik n x + B( ,-ik nX 

Ce ik n X + De -ik n X 



(88) 
(89) 



Applying the boundary conditions: ([86]) and ([87|l . one 
obtains: 

= A + B + C + D, (90) 

= Ae l{kn+K)L + De^ l{kn+K)L 

+ Be-' l{k "- K)L + Ce i( - k "- K)L . (91) 

The boundary conditions are satisfied with the choice: 

A=-D, B = C = 0, (92) 

which leads to sin[(fc n + K)L] = 0. Therefore the allowed 
values of k n are given by 



nil Air 

kn = - z 7, j 
L 6ao 

and the eigenenergies are given by: 



£ 2 — k y + k n 



No surface states exist in this case. 



(93) 



(94) 
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I. Dirac fermions in a magnetic field 

Let us now consider the problem of a uniform magnetic 
field B applied perpendicular to the graphene plane 2 . 
We use the Landau gauge: A = B(—y,0). Notice that 
the magnetic field introduces a new length scale in the 
problem: 



(95) 



which is the magnetic length. The only other scale in 
the problem is the Fermi-Dirac velocity. Dimensional 
analysis shows that the only quantity with dimensions of 
energy we can make is Vf/Ib- In fact, this determines 
the cyclotron frequency of the Dirac fermions: 



V2 V f 



(96) 



(the V2 factor comes from the quantization of the prob- 
lem, see below). Eqs. |96|) and |95|) show that the cy- 
clotron energy scales like y/B, in clear contrast with the 
non-relativistic problem where the cyclotron energy is 
linear in B. This implies that the energy scale associated 
with the Dirac fermions is rather different from the one 
find in the ordinary 2D electron gas. For instance, for 
fields of the order B w 10 T the cyclotron energy in the 
2D electron gas is of the order of 10 K. In contrast, for the 
Dirac fermion problem, for the same fields, the cyclotron 
energy is of the order of 1,000 K, that is, two orders of 
magnitude bigger. This has strong implications for the 
observ ation of the quantum H all effect at room temper- 
ature i Novoselov et all l2007h . Furthermore, for B = 10 
T the Zeeman energy is relatively small, g^sB w 5 K, 
and can be disregarded. 

Let us now consider the Dirac equation in more detail. 
Using the minimal coupling in (fl9|) (i.e., replacing — iV 
by — iV + eA/c) we find: 

v F [<?• (-iV + eA/c)] V(r) = E ip(r) , (97) 

in the Landau gauge the generic solution for the wave- 



function has the form ip(x,y) 
equation reads: 



e <t>{y), and the Dirac 





-dy-k + Bey/c 



dy — k+Bey/c 







4>{y) = Ed>{y\98) 



that can be rewritten as: 



O 
0t 



(99) 



! The problem of transver se magnetic and electric fields can also 
be so lved exactly. See: llLukose et all . 120071 : iPeres and Castrol . 
12007ft . 




Figure 18 (Color online) SdH oscillations observed in longi- 
tudinal resistivity p xx of graphene as a function of the charge 
carrier concentration n. Each peak corresponds to the pop- 
ulation of one Landau level. Note that the sequence is not 
interrupted when passing through the Dirac point, between 
electrons and holes. The period of oscillations An = 4B/<&o, 
where B is the applied field and $o is the flux quantum 
l|Novoselov et ad 120053 ). 



or equivalently: 

(Oa + 



= (2E/u c )<f>, 



(100) 



where a ± = cr x ± icr y , and we have defined the dimen- 
sionless length scale: 



i B k, 



and ID harmonic oscillator operators: 



(101) 



o 

C9t 



1 



(102) 



that obey canonical commutation relations: [0,0^] = 1. 
The number operator is simply: N = O^O. 

Firstly, we notice that (j 100p allows for a solution with 
zero energy: 



(o<7 + + oV)^ = 0, 



(103) 



and since the Hilbert space generated by a is of dimension 
2, and the spectrum generated by is bounded from 
below, we just need to ensure that: 



a~(f> 



(104) 



in order for (|103fl to be fulfilled. The obvious zero mode 
solution is: 



0o(O=Vo(O®l^>> 



(105) 
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where | JJ.) indicates the state localized on sublattice A 
and | ft) indicates the state localized on sublattice B. 
Furthermore, 



(106) 



is the ground states of the ID harmonic oscillator. All 
the solutions can now be constructed from the zero mode: 



&v,±(e) = ^-i(o®ift)±^(o 



(8 I V) 



and their energy is given by ijMcClurd . Il956f ) : 

E±(N) = ±u c VN, 



(107) 



(108) 



where N = 0,1,2,... is a positive integer, iPn{£,) is 
the solution of the ID Harmonic oscillator (explicitly: 
Vjv(0 = 2- N / 2 (M)- 1 / 2 exp{-e/2}H N (0 where H N (0 
is a Hermite polynomial). The Landau levels at the 
opposite Dirac point, K', have exactly the same spec- 
trum and hence each Landau level is doubly degener- 
ate. Of particular importance for the Dirac problem dis- 
cussed here is the existence of a zero energy state N = 
which is responsible, as we are going to show, to the 
anomalies observed in the quantum Hall effect. This 
particular Landau level structure has been observed by 
many different experimental probes, from Shubnikov-de 
Haas oscillations in single layer graphene ( see Fig. 1X8]) 



(Novoselov et al, 2005a; Z 



spectroscopy (IJiang et al 



lang et all |2005j) , to infrared 
2007alb and to scanning tun- 



neling spectroscopy ( Li and Andreil . I2007TJ (STS) on a 
graphite surface. 

J. The anomalous integer quantum Hall effect 

In the presence of disorde r Landau levels g et broadened 



in tne presence ot disorde r Landau levels g et broadened 
and mobility edges appear 1 Laughlinl Il98lh I . Notice that 
there will be a Landau level at zero energy that separates 
states with hole character (/i < 0) from states with elec- 
tron character (/x > 0). The components of the resistivity 
and conductivity tensors are related by: 



Pxx 
Pxy 



a 2 + a 2 

^ xx 1 xy 



a 2 + a 2 

xx 1 xy 



(109) 



where a xx (p xx ) is the longitudinal component and a xy 
(Pxy) is the Hall component of the conductivity (resistiv- 
ity). When the chemical potential is inside of a region 
of localized states the longitudinal conductivity vanishes, 
a xx = 0, and hence: p xx = 0, p xy = l/a xy . On the other 
hand, when the chemical potential is a region of delo- 
calized states, when the chemical potential is crossing a 
Landa u level, we have a x x ^ and u xy varies continu- 
ously llSheng et all I2006L 120071 ) . 




Figure 19 (Color online) Geometry of Laughlin's thought ex- 
periment with a graphene ribbon: a magnetic field B is ap- 
plied normal to the surface of the ribbon, a current / circles 
the loop, generating a Hall voltage Vh, and a magnetic flux 



The value of a xy in the region of localized states can 
be obtained from Laughlin's gauge invariance argument 
()Laughlinl . 1 198lh : one imagines making a graphene rib- 
bon such as the one in Fig. [19] with a magnetic field B 
normal through its surface and a current / circling its 
loop. Due to the Lorentz force the magnetic field pro- 
duces a Hall voltage Vh perpendicular to the field and 
current. The circulating current generates a magnetic 
flux $ that threads the loop. The current is given by: 



T 5E 



(110) 



where E is the total energy of the system. The localized 
states do not respond to changes in only the delocal- 
ized ones. When the flux is changed by a flux quantum 
(5$ = $ — hc/e the extended states remain the same 
by gauge invariance. If the chemical potential is in the 
region of localized states, all the extended states below 
the chemical potential will be filled both before and after 
the change of flux by $o- However, during the change of 
flux an integer number of states enter the cylinder at one 
edge and leave at the opposite edge. 

The question is how many occupied states are trans- 
ferred between edges. Let us consider a naive and as 
shown further incorrect calculation in order to show the 
importance of the zero mode in this problem. Each Lan- 
dau level contributes with one state times its degeneracy 
g. In the case of graphene we have g = 4 since there 
are 2 spin states and 2 Dirac cones. Hence, we would 
expect that when the flux changes by one flux quantum 
the change in energy would be SE lnc , = ±4NeVn, where 
N is an integer. The plus sign applies to electron states 
(charge +e) and the minus sign to hole states (charge 
— e). Hence, we would conclude that Ii nc . = ±4(e 2 / K)Vh 



21 



p„ (kQ) 



(4e 2 /h) 




is, at the Dirac point, it would predict a Hall plateau at 
N = with <J xy ,inc- = which is not possible since there 
is a N = Landau level, with extended states at this en- 
ergy. The solution for this paradox is rather simple: be- 
cause of the presence of the zero mode which is shared by 
the two Dirac points, there are exactly 2 x (2iV+ 1) occu- 
pied states that are transferred from one edge to another. 
Hence, the change in energy is SE — ±2(2iV + l)eVH for 
a change of flux of (5$ = hc/e. Therefore , the Ha l l con- 
ductivity is ( Gusvnin and Sharapoyl 120051 : iHerbutl . 120071 : 
IPeres et all 12006 J dl: ISchakell . Il99lh : 



I c SE e 2 

^ = 77- = 7^= ±2 ( 27V+1 )l7 
Vh Vh 5$ h 



(111) 



without any Hall plateau at N = 0. This amazing re- 
sult has been observed ex perimentally ( Novoselov et all 
l2005at IZhang et all l2005f ) as shown in Fig[20l 



Figure 20 (Color online) Quantum Hall effect in graphene as 
a function of charge carrier concentration. The peak at n = 
shows that in high magnetic fields there appears a Landau 
level at zero energy where no states exist in zero field. The 
field draws electronic states for this level from both conduc- 
tion and valence bands. The dashed line in dicate plateaus in 
a X y de scribed by Eq. (JTTTJ). Adapted from ijNovoselov et all 
l2005al ). 

and hence <J xy ,inc. = I /Vh = ±4Ne 2 /h, which is the 
naive expectation. The problem with this result is that 
when the chemical potential is exactly at half-filling, that 



K. Tight-binding model in a magnetic field 

In the tight-binding approximation the hopping inte- 
grals are replaced by a Peierls substitution: 

e ie &' A - dr i R ,R' = e*% A d V,R' . (H2) 
where £r,r' represents the hopping integral between the 
sites R and R', with no field present. The tight-binding 
Hamiltonian for a single graphene layer, in a constant 
magnetic field perpendicular to the plane, is conveniently 
written as, 



[ e I7r * ™ 2 a] J {m^n)b a (m 1 n)+e ™®a n al r (m,n)b a -(m—l,n~(l—z)/2)+e''™ <s> o n 2 aj.(m, n)&o-(m, n— z)+h.a], 

013) 



holding for a graphene stripe with a zigzag (z = 1) and Fourier transforming along the x direction gives, 
armchair (z = —1) edges oriented along the a;— direction. 



H = -tJ2 [e ftr ^ ni * £ ot(fc,n)6 (T (jfe ) n) + e~"^V fca 4 (k, n)b a (k, n - (1 - z)/2) + e"4™ £ r i a t (fc, »)& a (fc, n - z) + h.c.]. 



Let us now consider the case of zigzag edges. The 
eigenp roblem can be rewritten in terms of Harper's equa- 
tions (lHarnerl . 119551 ). and for zigzag edges we obtain 



(|Rammall . [l985h : 



£^ ife a(fc,n) 



-t\ e lka / 2 2cos(ir— n- — )/3(k,n) 



+ (3(k,n-l)] 
= -t[ e - lka/2 2 
+ a(k,n + 1)], 



-t[e- ifc °/ 2 2cos(7r— n - —)a(k 
1 v $ 2 ' y 



(114) 

n) 

(115) 
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zigzag: N=200, <K* =1/701 armchair. N=200, ^=1/701 




momentum ka momentum ka 



E(cV) E(cV) 




Figure 22 (Color online) Landau levels of the graphene stacks 
shown in Fig ll3l The applied magnetic field is 1 T. 



Figure 21 (Color online) Fourteen Energy levels of tight- 
binding electrons in graphene in the presence of a magnetic 
flux $ = $o/701, for a finite stripe with N = 200 unit cells. 
The bottom panels are zoom in images of the top ones. The 
dashed line represents the chemical potential jx. 



where the coefficients a(k,n) and f3(k,n) show up in 
Hamiltonian's eigenfunction, \tp(k)) } written in terms of 
lattice-position-state states as: 



\m) 



n,<j 



(a(k, n)\a; k, n, a) + f3(k, n)\b; A;, n, a}) 



(116) 

Eqs. I|114j) and (|1 15j) hold in the bulk. Considering that 
the zigzag ribbon has TV unit cells along its width, from 
n = Oton = N — 1, the boundary conditions at the edges 
are obtained from Eqs. (|1 14f) and (|115p . and read 

E^ k a{k, 0) = -te lka ' 2 2 cos Hjj (3(k, 0) , (117) 



potential. This counting produces the value (2N + 1), 
with N = 0,1,2,... (for the case of Fig. [21] one has 
N = 2). From this counting the Hall conductivity is 
given, including an extra factor of two accounting for the 
spin degree of freedom, by 



±2— (2N- 



1) = ±4- 



N 



(119) 



E^ k /3(k, N-l)=-2te- lka/2 cos 



' ka 
$o 2 



Eq. I|119p leads to the anomalous integer quantum Hall 
effect discussed previously, which is the hallmark of single 
layer graphene. 



L. Landau levels in graphene stacks 

The dependence of the Landau level energies on the 
Landau index TV roughly follows the dispersion relation 
of the bands, as shown in Fig. [22j Note that, in a trilayer 
with Bernal stacking, two sets of levels have a s/N depen- 
dence, while the energies of other two depend linearly on 
N. In an infinite rhombohedral stack, the Landau lev- 



remain discrete and quasi-2D ijGuinea et all 1200 



Similar equations hold for a graphene ribbon with arm- 
chair edges. 

In Fig. [21] we show fourteen energy levels, around zero 
energy, for both the zigzag and the armchair cases. The 
formation of the Landau levels is signaled by the presence 
of flat energy bands, following the bulk energy spectrum. 
From Fig. [21] it is straightforward to obtain the value of 
the Hall conductivity in the quantum Hall effect regime. 
Let us assume that the chemical potential is in between 
two Landau levels at positive energies, as represented by 
the dashed line in Fig. [21] The Landau level structure 
shows two zero energy modes, one of them is electron-like 
(hole-like) , since close to the edge of the sample its energy 
is shifted upwards (downwards) . The other Landau levels 
are doubly degenerate. The determination of the values 
for the Hall conductivity is done by counting how many 
energy levels (of electron-like nature) are below chemical 



Note that, even in an infinite stack with the Bernal struc- 
ture, the Landau level closest to E = forms a band 
which does not overlap with the other Landau levels, 
leadi ng to the possib i lity o f a 3D integer quantum Hall ef- 
fect (iBernevig et aZ.l.l2007l:lKopelevich et aLl . [200^. l2003t 
iLuk'vanchuk and Kopelevichl . 2004 ). 

The optical transitions between Landau levels can also 
be calculated. The selection rules are the same as for a 
graphene single layer, and only transitions between sub- 
bands with Landau level indices M and N such that 
\N\ = |A/±1| are allowed. The resulting transitions, with 
their respective spectral strengths, are shown in Fig. [23] 
The transitions are grouped into subbands, which give 
rise to a continuum when the number of layers tends to 
infinity. In Bernal stacks with an odd number of lay- 
ers, the transitions associated to Dirac subbands with 
linear dispersion have the largest spectral strength, and 
they give a significant contribution to the total absorp- 
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Figure 23 (Color online) Relative spectral strength of the low 
energy optical transitions between Landau levels in graphene 
stacks with Bernal ordering and an odd number of layers. The 
applied magnetic field is 1 T. Top left: 3 layers. Top right: 
11 layers. Bottom: 51 layers. The large red circles are the 
transitions in a single layer. 



tion even if the numbe r of layers is large, Ni, < 30 
()Sadowski et aZ.I . l2006h . 
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M. Diamagnetism 



Back in 1952 Mrozowski (jMrozowskil . llQ52h stud- 
ied diamagnetism of polycrystalline graphite and other 
condensed-matter molecular-ring systems. It was con- 
cluded that in such ring systems diamagnetism has two 
contributions: (1) diamagnetism from the filled bands of 
electrons; (2) Landau diamagnetism of free electrons and 
holes. For graphite the second source of diamagnetism is 
by far the l argest of the two . 

McClure ( McClurel . Il956h computed diamagnetism of 
a 2D ho neycomb l attice using the theory introduced by 
Wallace l| Walked . fl947h . a calcula tion he la t er ge neral- 
ized to three dimensional graphite ( McClure] . Il960f ). For 
the honeycomb plane the magnetic susceptibility in units 
of emu/g is 



X 



0.0014 



7g sech 2 



/' 



2k B T 



(120) 



where /i is the Fermi energy, T the temperature, and 
ks the Boltzmann constant. For graphite the magnetic 
susceptibility is anisotropic and the difference between 
the susceptibility parallel to the principal axis and that 
perpendicular to the principal axis is -21. 5x 10 -6 emu/g. 
The susceptibility perpendicular to the principal axis is 
about the free-atom susceptibility of -0.5xl0 -6 emu/g. 
In the 2D model the susceptibility turns out to be large 
due to the presence of fast moving electrons, with a ve- 
locity of the order of vf — 10 6 m/s, which in turn is a 
consequence of the large value of the hopping parameter 
7o . In fact the susceptibility turns out to be proportional 
to the square of 70 . Later Sharma et al. extended the cal- 
culation of McClure for graphite by including the effect 



of trigonal warping and s howed that the low t emperature 
diamagnetism increases jlSharma et al. . 1973) . 

Safran and DiSalvo ijSafran and DiSalvol . I19791L in- 
terested in the susceptibility of graphite intercalation 
compounds, recalculated, in a tight-binding model, the 
susceptibility perpe ndicular to a gr aphite plane using 
Fukuyama's theory 1 Fukuy amal . 1 1 9 7 lh . which includes in- 
terband transitions . The results a gree with those pub- 
lished by McClure (|McClurdJl956i ). Later, Safran com- 
puted the susceptibility of a graphene bilayer showing 
that this system is diamagnetic at small values of the 
Fermi energy, but there appears a paramagnetic peak 
when the Fermi energy is of the order of the interlayer 
coupling (|Safranl . [l98l . 

The magnetic susceptibility of other carbon based ma- 
terials, as carb on nanotubes and Cap molecular solids 
was measured ( Heremans et all . Il994l ) showing a dia- 
magnetic response at finite magnetic fields different from 
that of graphite. The study of the magnetic response 
of nanographite ribbons with both zig-zag and arm-chair 
edges was done by Wakabayashi et al. using a numeri- 



cal d ifferentiation of the free energy (jWakabayashi et ali . 
1999). From these two systems, the zig-zag edge state 
is of particular interest since the system supports edge 
states even in the presence of a magnetic field, leading to 
very high density of states near the edge of the ribbon. 
The high temperature response of these nanoribbons was 
found to be diamagnetic whereas the low temperature 
susceptibility is paramagnetic. 

The Dirac-like nat ure of the el e ctroni c quasiparti- 
cles in graphene led ( Ghosal et all . l2007h to consider 
in general the problem of the diamagnetism of nodal 
fermions and Nakamura to study the orbital magnetism 
of D irac fermions in weak magnetic fields 1 Nakamural . 
l2007l L Koshino and Ando considered the diamagnetism 
of disordered graphene in the self consistent Born approx- 
imation, wit h a disorder potential of the form U{r) = 
lui5(r - R) (jKoshino and Andol . hoOlh . At the neutral- 
ity point and zero temperature the susceptibility of dis- 
ordered graphene is given by 



X(0) 



2W 

s^^-iV 



(121) 



where g s = g v — 2 is the spin and valley degeneracies, W 
is a dimensionless parameter for the disorder strength, 
defined as W = niuf / 'Attjq , rij the impurity concentra- 
tion, and To is given by Tq = e c exp[— 1 / (2W)] with e c 
a parameter defining a cut-off function used in the the- 
ory. At finite Fermi energy ep and zero temperature the 
magnetic susceptibility is given by 



9v9s 



2 2 w 

3tt 6 70 \e F \ 



(122) 



In t he clean limit t he sus cep tibility i s given 
bv dKoshino and Andol . l2007t iMcClurd . Il956t 
ISafran and DiSalvcLll97gl : 



x(cf) 



9v9s 
6ir 



(123) 
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N. Spin orbit coupling 

Spin-orbit coupling describes a process in which an 
electron changes simultaneously its spin and angular mo- 
mentum or, in general, moves from one orbital wavefunc- 
tion to another. The mixing of the spin and the orbital 
motion is a relativistic effect, which can be derived from 
Dirac's model of the electron. It is large in heavy ions, 
where the average velocity of the electrons is higher. Car- 
bon is a light atom, and the spin orbit interaction is ex- 
pected to be weak. 

The symmetries of the spin orbit interaction, however, 
allow the formation of a gap at the Dirac points in clean 
graphene. The spin orbit interaction leads to a spin de- 
pendent shift of the orbitals, which is of a different sign 
for the two sublattic es, acting as an effective mass within 
each Dirac po i nt (iDresselhaus and Dresselhausl. 1 19651: 
iKane and Meld . 120051 ; I Wang and ChakrabortyL l2007al ). 
The appearance of this gap leads to a non trivial spin 
Hall conductance, in similar way to other models which 
study the parity a nomaly in relati vistic field theory in 
(2+1) dimensions ( Haldanel . Il988h . When the inver- 
sion symmetry of the honeycomb lattice is broken, ei- 
ther because the graphene layer is curved or because 
an external electric field is applied (Rashba interaction) 
additional terms, which m odulat e the nearest neigh- 
bor hopping, are induced ( Andd . l2000l ). The intrin- 
sic and extrinsic spin orbit interactions can be written 
as (IDresselhaus and Dresselhausl . Il965l : IKane and Meld . 
120051 ): 



Hsount = A so J d 2 r¥ {r)s z a z f z * (r) , 

H SO ;e Xt = \r Jd 2 r¥ (r)(-s x a y +s y a x T z )i>(r) ,(124) 



where a and f are Pauli matrices which describe the sub- 
lattice and valley degrees of freedom, and s are Pauli ma- 
trices acting on actual spin space. A so is the spin-orbit 
coupling and Xr is the Rashba coupling. 

The magnitude of the spin orbit coupling in graphene 
can be inferred from the known spin orbit coupling in 
the carbon atom. This coupling allows for transitions 
between the p z and p x and p y orbitals. An electric field 
induces also transitions between the p z and s orbitals. 
These intra-atomic processes mix the ir and a bands 
in graphene. Using second order perturbation theory, 
one obtains a coupling betwe en the low energy states 
in the 7r band. Tight-binding (jHuertas-Hernando et all 
l2006l ; IZarea and Sandler| . l2007l) and band structure calcu- 
lations ijMin et art l2006l : lYao et aZ.I . [2007h give estimates 
for the intrinsic and extrinsic spin-orbit interactions in 
the range 0.01 — 0.2 K, and hence much smaller than the 
other energy scales in the problem (kinetic, interaction, 
and disorder). 



III. FLEXURAL PHONONS, ELASTICITY, AND 
CRUMPLING 

Graphite, in the Bernal stacking configuration, is a lay- 
ered crystalline solid with 4 atoms per unit cell. Its basic 
structure is essentially a repetition of the bilayer struc- 
ture discussed earlier. The coupling between the lay- 
ers, as we discussed, is weak and, therefore, graphene 
has been always the starting point for the d iscussion of 
phonons in graphite i Wirtz and Rubiol . l2004h . Graphene 
has two atoms per unit cell and if we consider graphene as 
strictly 2D it should have 2 acoustic modes (with disper- 
sion uj ac (k) oc k as k — > 0) and 2 optical modes (with dis- 
persion w op (A;) oc constant, as k — > 0) solely due to the in- 
plane translation and stretching of the graphene lattice. 
Nevertheless, graphene exists in the 3D space and hence 
the atoms can oscillate out-of-plane leading to 2 out-of- 
plane phonons (one acoustic and another optical) called 
flexural modes. The acoustic flexural mode has disper- 
sion Wfl ox (fc) oc k 2 as k — > which represents the transla- 
tion of the whole graphene plane (essentially a one atom 
thick membrane) in the perpendicular direction (free par- 
ticle motion). The optical flexural mode represents the 
out-of-phase out-of-plane oscillation of the neighboring 
atoms. In first approximation, if we neglect the cou- 
pling between graphene planes, graphite has essentially 
the same phonon modes, albeit they are degenerate. The 
coupling between planes has two main effects: (1) it lifts 
the degeneracy of the phonon modes, and (2) leads to a 
strong suppression of the energy of the flexural modes. 
Theoretically, flexural modes should become ordinary 
acoustic and optical modes in a fully covalent 3D solid, 
but in practice, the flexural modes survive due to the 
fact the planes are coupled by weak van der Waals-like 
forces. Thes e modes have been mea sured experimentally 
in graphite ( Wirtz and Rubiol . l2004h . Graphene can also 
be obtained as a suspended membrane, that is, without 
a substrate, being supp orted only by a scaffold or bridg- 



ing m icron-scale gaps l|Bunch et all l2007t iMever et all 
l2007allbh . Figure [24] shows a suspended graphene sheet 
and an atomic resolution image of its crystal lattice. 

Because the flexural modes disperse like k 2 they dom- 
inate the behavior of structural fluctuations in graphene 
at low energies (low temperatures) and long wave- 
lengths. It is instructive to understand how these 
mod es appear from the poi n t of view of elastic i ty the - 
ory ijChaikin and Lubenskvl . Il995t iNelson et all . l2004h . 
Consider, for instance, a graphene sheet in 3D and let 
us parameterize the position of the sheet relative of a 
fixed coordinate frame in terms of the in-plane vector r 
and the height variable h(r) so that a position in the 
graphene is given by the vector R = (r, h(r)). The unit 
vector normal to the surface is given by: 



N = 



- Vh 



(125) 



where V = (d x ,d y ) is the 2D gradient operator, and 
z is the unit vector in the third direction. In a flat 
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to the elastic energy: 



Figure 24 (Color online) Suspended graphene sheet, (a) 
Bright-field transmission-electron-microscope image of a 
graphene membrane. Its central part (homogeneous and 
featureless region) is monolayer graphene. Adapted from 
ijMever et a/.l . l2007ah . (b) Despite being only one atom thick, 
graphene remains a perfect crystal as this atomic resolution 
image shows. The image is obtained in a scanning transmis- 
sion electron microscope. The visible peri odicity is given b y 
the lattice of benzene rings. Adapted from l|Nair et a/.Ll2008l ). 



graphene configuration all the normal vectors are aligned 
and therefore V • N = 0. Deviations from the flat config- 
uration requires misalignment of the normal vectors and 
costs elastic energy. This elastic energy can be written 
as: 



Eo = - d 2 r (V • N) 



d 2 r(V 2 /i) 2 (126) 



where k is the bending stiffness of graphene, and the 
expression in terms of h(r) is valid for smooth distortions 
of the graphene sheet ((V/i) 2 <C 1). The energy (j 1 26[) is 
valid in the absence of a surface tension or a substrate 
which break the rotational and translational symmetry 
of the problem, respectively. In the presence of tension 
there is an energy cost for solid rotations of the graphene 
sheet (V/i 7^ 0) and hence a new term has to be added 



E T = 1 d 2 r(Whf 



(127) 



where 7 is the interfacial stiffness. A substrate, described 
by a height variable s(r), pins the graphene sheet through 
van der Waals and other electrostatic potentials so that 
the graphene configuration tries to follow the substrate 
h(r) ~ s(r). Deviations from this configuration cost ex- 
tra elastic energ y that can be approx i mated by a har- 
monic potential 1 Swain and Andelmanl . [l999h : 



Ea = l j ' d 2 v{s{v)-h(v)) 2 , 



(128) 



where v characterizes the strength of the interaction po- 
tential between substrate and graphene. 

Firstly, let us consider the free floating graphene prob- 
lem l|126|) that we can rewrite in momentum space as: 



(129) 



We now canonically quantize the problem by introducing 
a momentum operator Pk that has the following commu- 
tator with /ik: 



[hk,Pk>] = i$k,k 
and write the Hamiltonian as: 

£{- kPk 

k 



H 



2a 



— h-khk 



(130) 



(131) 



where a is graphene's 2D mass density. From the Heisen- 
berg equations of motion for the operators it is trivial 
to find that oscillates harmonically with a frequency 
given by: 



Wflc X (k) = ( - 



1/2 



(132) 



which is the long wavelength dispersion of the flexural 
modes. In the presence of tension it is easy to see that 
the dispersion is modified to: 



w (k) = kj-k 2 + - 



(133) 



indicating that the dispersion of the flexural modes be- 
comes linear in k, as k — > 0, under tension. That is what 
happens in graphite where the interaction between layers 
breaks the rotational symmetry of the graphene layers. 

Eq. (|132p also allows us to relate the bending energy 
of graphene with the Young modulus, Y, of graphite. 
The fundamental resonance frequency of a macroscopic 
graph ite sample of thickness t is given by 1 Bunch et all 
l2r)07f i: 



'7 



1/2 



tk 2 



(134) 
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where p = ajt is the 3D mass density. Assuming that 
(|134|) works down to the single plane level, that is, when 
t is the distance between planes, we find: 



(135) 



which provides a simple relationship between the bending 
stiffness and the Young modulus. Given that Y w 10 12 
N/m and t m 3.4 A we find, k w 1 eV. This result is in 
good agr eement with ab i nitio calculations of the bending 
rigidity ijLenoskv et all Il992l : iTu and Ou-Yand. l2002fl 
and experiments in graphene resonators IjBunch et all 
l2007h . 

The problem of structural order of a "free floating" 
graphene sheet can be fully understood from the exis- 
tence of the flexural modes. Consider, for instance, the 
number of flexural modes per unit of area at a certain 
temperature T: 



N, 



d 2 k l_ 

(2^) 2 ™ k ~ 2ir 



dk- 



1 



(136) 



where is the Bose-Einstein occupation number (/? = 
l/(fcsT)). For T^O the above integral is logarithmically 
divergent in the infrared (k — > 0) indicating a divergent 
number of phonons in the thermodynamic limit. For a 
system with finite size L the smallest possible wave vector 
is of the order of fc min ~ 2n/L. Using fc m ; n as a lower cut- 
off in the integral (|136[) we find: 



JVph = 5: ln (i- e -^ 



where 



27T /K\!/4 



(137) 



(138) 



is the thermal wavelength of the flexural modes. Notice 
that that when L ^> Lt the number of flexural phonons 
in l|137jl diverges logarithmically with the size of the sys- 
tem: 



2tt 



L 



(139) 



indicating that the system cannot be structurally or- 
dered at any finite temperature. This is noth- 
i ng but the crumpling instability of soft memb ranes 
(jChaikin and LubenskyL Il995t iNelson et all l2004 l- For 
L <C Lt one finds that N p h goes to zero exponentially 
with the size of the system indicating that systems with 
finite size can be flat at sufficiently low temperatures. 
Notice that for k w 1 eV, p w 2200 kg/m 3 , t = 3.4 
A(er « 7.5 x 10~ 7 kg/m 2 ), and T w 300 K, we find 
Lt »li indicating that free floating graphene should 
always crumple at room temperature due to thermal 
fluctuations associated with flexural phonons. Never- 
theless, the previous discussion only involves the har- 
monic (quadratic part) of the problem. Non-linear effects 



such as large bending deformations i Peliti and LeibleiL 
Il985h . the coupling between flexural and in-plane modes 
(or phonon-phonon interact i ons I Bonini et all l2007t 
iRadzihovskv and Le Doussall. 1992?) ) and the presence 
of topological defects ( Nelson and Peliti Il987h can lead 
to strong renormalizations of the bending rigidity, driv- 
ing t he system toward a flat pha se at low tempera- 
tures (jChaikin and LubenskyL Il995l ) . This situation has 
been c onfirmed in numerical s i mulations of free grap hene 
sheets ijAdebpour et all I20071 : iFasolino et all 120071 1 . 

The situation is rather different if the system is un- 
der tension or in the presence of a substrate or scaffold 
that can hold the graphene sheet. In fact, static rippling 
of graphene flakes suspended on scaffolds have been ob- 
served for single layer as well as bilayers 1 Meyer et all 
l2007allhh . In this case the dispersion, in accordance with 
(|133p . is at least linear in k, and the integral in (|136|) 
converges in the infrared (k — > 0) indicating that the 
number of flexural phonons is finite and graphene does 
not crumple. We should notice that these thermal fluctu- 
ations are dynamic and hence average to zero over time, 
therefore, the graphene sheet is expected to be flat under 
these circumstances. Obviously, in the presence of a sub- 
strate or scaffold described by (j 1 28fl static deformations 
of the graphene sheet are allowed. Also, hydrocarbon 
molecules that are often present on top of free hanging 
graphene membranes could quench flexural fluctuations 
making them static. 

Finally, one should notice that in the presence of a 
metallic gate the electron-electron interactions lead to 
the coupling of the phonon modes to the electronic ex- 
citations in the gate. This coupling could be partially 
responsible to the damping of the phono n modes due to 
dissipative effects ( Seoanez et all l2007Ll as observed in 
graphene resonators ijBunch et all 2007t l. 



IV. DISORDER IN GRAPHENE 

Graphene is a remarkable material from the electronic 
point of view. Because of the robustness and specificity 
of the sigma bonding, it is very hard for alien atoms to 
replace the carbon atoms in the honeycomb lattice. This 
is one of the reasons why the electron mean free path 
in graphene can be so long, reaching up to one microm- 
eter in the existing samples. Nevertheless, graphene is 
not immune to disorder and its electronic properties are 
controlled by extrinsic as well as intrinsic effects that are 
unique to this system. Among the intrinsic sources of 
disorder we can highlight: surface ripples and topologi- 
cal defects. Extrinsic disorder can come about in many 
different forms: adatoms, vacancies, charges on top of 
graphene or in the substrate, and extended defects such 
as cracks and edges. 

It is easy to see that from the point of view of single 
electron physics (that is, terms that can be added to (J3J), 
there are two main terms that disorder couples to. The 
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first one is a local change in the single site energy, 

H dd = J2 V > ( a l a * + b t b i) ' ( 14 °) 

i 

where Vi is the strength of the disorder potential on site 
Ri, which is diagonal in the sublattice indices and hence, 
from the point of view of the Dirac Hamiltonian lfl8|) , can 
be written as: 

H dd = f d 2 r ]T 7 (r)*t(r)* a (r), (141) 

■* 0=1,2 

which acts as a chemical potential shift for the Dirac 
fermions, that is, shifts locally the Dirac point. 

Because of the vanishing of the density of states 
in single layer graphene, and by consequence the lack 
of electrostatic screening, charge potentials may be 
rather important in determining the spectroscopic and 
transport properties (lAdam et all . 120071 : lAndol . l2006ri 
INomura and MacDonaldl . 12007 ). Of particular impor- 
tance is the Coulomb impurity problem where, 



V a (r) 



eo r 



(142) 



where eo is the dielectric constant of the medium. The 
solution of the Dirac equation for the Co ulomb poten- 
tial in 2D can be studied a nalytic a lly dBiswas et al 



2007: DiVincenzo an 



Pereira et all . l2007bl : ISlivtov et all . \200H) . Its solution 



dMeld. Il98l iNovikovl 12007a: 



has many of the f eatures of th e 3D relativistic hydro- 
gen atom problem i Bavml . fl969). Just as in the case of 
the 3D problem the nature of the eigenfunctions depends 
strongly on graphene's dimensionless coupling constant: 



9 



Ze 2 
e v F 



(143) 



Notice, therefore, that the coupling constant can be var- 
ied by either changing the charge of the impurity, Z, 
or modifying the dielectric environment and changing 
eo. For g < g c = 1/2 the solutions of this problem are 
given in terms of Coulomb wavefunctions with logarith- 
mic phase shifts. The local density of states (LDOS) is 
affected close to the impurity due the electron-hole asym- 
metry generated by the Coulomb potential. The local 
charge density decays like 1/r 3 plus fast oscillations of 
the order of the lattice spacing (in the continuum limit 
this would give rise t o a Dir ac delta function for the den- 
sity (jKolezhuk et all |2006( )). Just like in 3D QED, the 
2D problem becomes unstable for g > g c = 1/2 leading 
to super-critical be havior and the so-ca l led fa ll of elec- 
tron to the center rtLandau and Lifshitd . fl98lh . In this 
case the LDOS is strongly affected by the presence of the 
Coulomb impurity with the appearance of bound states 
outsid e the band and scat tering resonances within the 
band I Pereira et al ]. l2007hl ) and the local electronic den- 
sity decays monotonically like 1/r 2 at large distances. 



It has been argued I Schedin et all . l2007t ) that with- 
out high vacuum environment these Coulomb effects can 
be strongly suppressed by large effective dielectric con- 
stants due to the p resence of a nanom eter thin layer 
of absorbed water (jSabio et all l2007t ). In fact, ex- 



periments in ultra-high vacuum conditions (|Chen et all . 
l2007bl ) display strong scattering features in the trans- 
port that can be associated to charge impurities. 
Screening effects that affect the strength and range 
of the Coulomb interaction, are rather non-trivial in 



graphene l|Fogler et all . l2007bl : IShklovskil 120071 ) and, 
theref ore, important fo r the in t erpretation of transport 
data dBardarson et all. l2007t iLewenkopf et all 120071 : 
INomura et all . 120071 : ISan- Jose et all . l2007[ ). 

Another type of disorder is the one that changes the 
distance or angles between the p z orbitals. In this case, 
the hopping energies between different sites are modified 
leading to a new term to the original Hamiltonian (0: 



id 

+ ^(afa, +&!&,)} , (144) 
or in Fourier space: 

H od = ^atfo k ^^ ab) e^ k - k ')- R -^- k '+h.c. 

+ (4a k <+&jA') E ^ aa) e^ k - k ')- R -^- k '(145) 



i.S a 



where St (St 



f is the change of the hopping en- 
ergy between orbitals on lattice sites R^ and Rj on the 
same (different) sublattices (we have written Rj = R; + 5 

where S a b is the nearest neighbor vector, and S aa is the 
next nearest neighbor vector). Following the procedure 
of Sec. III. Bl we project out the Fourier components of the 
operators close to the K and K' points of the BZ using 
(fT7"l) . If we assume that SUj is smooth over the lattice 
spacing scale, that is, it does not have an Fourier com- 
ponent with momentum K — K' (so the two Dirac cones 
are not coupled by disorder) , we can rewrite l|145|) in real 
space as: 



dfr „4(r)a](r)6i(r) + h.c. 



+ (f>(r) (a t 1 (r)a 1 (r) + & t 1 (r)& 1 (r))j , (146) 

with a similar expression for the cone 2 but with A re- 
placed by *4*, where, 



A(r) = J2 Sfiab) (r)t 
0(r) = J2 Stiaa Hr)< 



(147) 
(148) 



Notice that whereas 0(r) = 0*(r), because of the in- 
version symmetry of the two triangular sublattices that 



28 



make up the honeycomb lattice, A is complex because of 
lack of inversion symmetry for nearest neighbor hopping. 
Hence, 



A(r) =A x (r)+iA y (r). 



(149) 



In terms of the Dirac Hamiltonian l|T8[) we can rewrite 
(fT46l) as: 



d 2 r 



*i(r)cr • ^4(r)*i(r) 



(^(^(r) 



(150) 



where A — (A x ,A y ). This result shows that changes in 
the hopping amplitude lead to the appearance of vector, 
A, and scalar, <f>, potentials in the Dirac Hamiltonian. 
The presence of a vector potential in the problem indi- 
cates that an effective magnetic field B = c/(ev F )V x A 
should also be present, naively implying a broken time re- 
versal symmetry, although the original problem was time 
reversal invariant. This broken time reversal symmetry 
is not real since (|150|1 is the Hamiltonian around only 
one of the Dirac cones. The second Dirac cone is re- 
lated to the first by time reversal symmetry indicating 
that the effective magnetic field is reversed in the second 
cone. Therefore, there is no global broken symmetry but 
a compensation between the two cones. 



A. Ripples 

Graphene is a one atom thick system, the extreme 
case of a soft membrane. Hence, just like soft mem- 
branes, it is subject to distortions of its structure either 
due to thermal fluctuations (as we discussed in Sec. IIIII) 
or interaction with a substra te, scaffold, and absorbands 
( Swain and Andelmanl . 1 1 99*91 ) . In the first case the fluc- 
tuations are time dependent (although with time scales 
much longer than the electronic ones), while in the sec- 
ond case the distortions act as quenched disorder. In 
both cases, the disorder comes about because of the mod- 
ification of the distance and relative angle between the 
carbon atoms due to the bending of the graphene sheet. 
This type of off-diagonal disorder does not exist in ordi- 
nary 3D solids, or even in quasi-lD or quasi-2D systems, 
where atomic chains and atomic planes, respectively, are 
embedded in a 3D crystalline structure. In fact, graphene 
is also very different from other soft membranes because 
it is (semi) metallic, while previously studied membranes 
were insulators. 

The problem of the bending of graphitic systems and 
its effect on the hybridization of the tt orbitals has 
been studied a g reat deal in the con text of classical 
minimal surfaces (jLenoskv et all . 1992f) and app l ied to 
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fullerenes and carbon nanotubes (IKane and Mole, 
Tersofj. Il992t iTu and Ou-Yand . l2002t IXin et al 
Zhong-can et all . 1997wl ). In order to understand the 
effect of bending on graphene, consider the situation 



2000; 



shown in Figj25l The bending of the graphene sheet has 
three main effects: the decrease of the distance between 
carbon atoms, a rotation of the pz orbitals (compres- 
sion or dilation of the lattice are energetically costly due 
to the large spr ing constant of graphene « 57 eV/A 2 
(|Xin et aZ.l . l2000h1. and a re-hybridization betwee n tt and 
a orbitals l|Eun-Ah Kim and Castro Netol . l2007h . Bend- 
ing by a radius R decreases the distance between the 
orbitals from I to d = 2Rsm[£/(2R)} w I - l 3 /(24R 2 ) 
for R ^> L The decrease in the distance between 
the orbitals increas es the ove r lap b etween the two lobes 
of the p z orbital (|Harrisonl . Il980h : V ppa w V° pa [l + 



£ 2 /(12i? 2 )], where a = tt, a, and V p ° pa is the overlap for 
a flat graphene sheet. The rotation of the pz orbitals 
can be understood within the Slater-Koster formalism, 
namely, the rotation can be decomposed into a p z — p z 
(tt bond) plus a p x — p x (a bond) h ybridization wit h 
energies V pp7r and V ppa , respectively (jHarrisonl . ll98dT l : 



V(6) 



V m 



'-' 2 '^ - V ppa sin 2 (6) 



\0) 



V m 



V PP<T )(£/ (2i?)) 2 , leading to a decrease in the overlap. Fur- 
thermore, the rotation leads to re-hybridization between 
tt and a orbitals leading to a further shift in energ y 
of the order of ljEun-Ah Kim and Castro" Netol . 120071 1: 

fc* « (V 2 + V* )/(e n - e„). 




Figure 25 Bending of the surface of graphene by a radius R 
and its effect on the p z orbitals. 

In the presence of a substrate, as we discussed 
in SecHnl elasticity theory predicts that graphene 
can be expected to follow the substrate in a smooth 
way. Indeed, by minimizing the elastic energy (|126|1 . 
1271) . and (11281) with respe ct to the height h we get 



Swain and Andelmanl . fl~999T ) : 



KV 4 ft(r) - 7 V 2 /i(r) + vh(r) = vs(r) , (151) 
that can be solved by Fourier transform: 

*(k) 



h(k) 



i + {itkf + (e c k) 4 



(152) 



where 



7 \ 1/2 

vJ 

1/4 



(153) 
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Eq. I|152p gives the height configuration in terms of the 
substrate profile, and i t and £ c provide the length scales 
for elastic distortion of graphene on a substrate. Hence, 
disorder in the substrate translates into disorder in the 
graphene sheet (albeit restricted by elastic constraints). 
This picture ha s been confirm e d by STM measurements 
on g raphene (jlshigami et al\ . 120071 : IStolyarova et al\ . 
l2007t i in which strong correlations were found between 
the roughness of the substrate and the graphene topog- 
raphy. Ab initio band structure calculations als o give 
support to this scenario (jDharma- WardanaLl2007h . 

The connection between the ripples and the electronic 
problem comes from the relation between the height field 
h(r) and the local curvature of the graphene sheet, R: 



R(r) 



V 2 h{r) 



(154) 



and, hence we see that due to bending the electrons are 
subject to a pote ntial which depends on the structure o f 
a graphene sheet l|Eun-Ah Kim and Castro Netol . l2007t l : 



V(r) wl/°-a(V 2 /i(r)) 2 



(155) 



where a (a w 10 eV A 2 ) is the constant that depends on 
microscopic details. The conclusion from (|155j) is that the 
Dirac fermions are scattered by ripples of the graphene 
sheet through a potential which is proportional to the 
square of the local curvature. The coupling between ge- 
ometry and electron propagation is unique to graphene, 
and res ults in additional scattering a nd resistivity due to 
ripples (jKatsnelson and Geiml . l2008h . 



B. Topological lattice defects 

Structural defects of the honeycomb lattice like pen- 
tagons, heptagons and their combinations such as Stone- 
Wales defect (a combination of two pentagon-heptagon 
pairs) are also possible in graphene and can lead to scat- 
tering ( Cortijo and Vozmedianol . l2007af bh. These defects 
induce long range deformations, which modify the elec- 
tron trajectories. 

Let us consider first a disclination. This defect is equiv- 
alent to the deletion or inclusion of a wedge in the lattice. 
The simplest one in the honeycomb lattice is the absence 
of a 60° wedge. The resulting edges can be glued in such 
a way that all sites remain three-fold coordinated. The 
honeycomb lattice is recovered everywhere, except at the 
apex of the wedge, where a fivefold ring, a pentagon, is 
formed. One can imagine a situation where the nearest 
neighbor hoppings are unchanged. Nevertheless, the ex- 
istence of a pentagon implies that the two sublattices of 
the honeycomb structure can no longer be defined. A 
trajectory around the pentagon after a closed circuit has 
to change the sublattice index. 

The boundary conditions imposed at the edges of a 
disclination are sketched in Fig. [26l showing the identi- 
fication of sites from different sublattices. In addition, 




Figure 26 (Color online) Sketch of the boundary conditions 
associated to a disclination (pentagon) in the honeycomb lat- 
tice. 



the wavefunctions at the K and K' points are exchanged 
when moving around the pentagon. 

Far away from the defect, a slow rotation of the 
components of the spinorial wavefunction can be de- 
scribed by a gauge fi eld which acts on the valle y 
and sublattice indices ()Gonzalez et all . Il99l Il993bl) . 
This gauge field is technically non-abelian, although 
a transformation can be defined which makes the re- 
sulting Dirac equation e quivalent to one with a n effec- 
tive abelian gauge field (jGonzalez et~al\ . Il993bl ). The 
final continuum equation gives a reasonable descrip- 
tion of t he electronic spectrum of fu llerenes of differ- 
ent sizes ([Gonzalez et a^.l . ll99lfl993rl. and other struc- 
tures which contain pentagons (iKolesnikov and Osipov , 
2004 120061: [Lammert and Crespil . 12004 iLeClaid . I200C : 



Osipov et all 20031 ) . It is easy to see that an heptagon 



leads to the opposite effective field. 

An in-plane dislocation, that is, the inclusion of a semi- 
infinite row of sites, can be considered as induced by 
a pentagon and a heptagon together. The non-abelian 
field described above is canceled away from the core. A 
closed path implies a shift by one (or more) lattice spac- 
ings. The wavefunctions at the K and K 1 points acquire 



phases, e 



±2-rri/3 



, under a translation by one lattice unit. 



Hence, the description of a dislocation in the continuum 
limit requires an (abelian) vortex of charge ±27r/3 at its 
core. Dislocations separated over distances of the order 
of d lead to an ef fective flux through an area of perimeter 
I of the order of ( Morpurgo and Guinea! . l2006h : 



(156) 



where is the Fermi vector of the electrons. 

In general, a local rotation of the axes of the hon- 
eycomb lattice induces changes in the hopping which 
lead to mixing of the K and K' wavefunctions, lead- 
ing to a gauge field l ike th e one induced by a pen- 
tagon ( Gonzalez et all . l200lh . Graphene samples with 
disclination s and dislocati o ns ar e feasible in particular 
substrates ( Couraux et all . l2008f h and gauge fields re- 
lated to the local curvature are then expected to play 
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a crucial role in such structures. The resulting elec- 
tronic structure can be analyzed using the theory of quan- 
tum mechanics in curved space (Birrell and Daviesl 19821: 
ICortiio and Vozmediand . l2007al |bl; Ide Juan et all l2007t ). 



C. Impurity states 

Point defects, similar to impurities and vacancies, can 
nucleate electronic states in their vicinity. Hence, a con- 
centration of rii impurities per carbon atom leads to a 
change in the electronic density of the order of n, . The 
corresponding shift in the Fermi energy is ep ~ v-p^/ril. 
In addition, impurities lead to a finite elastic mean free 

— 1/2 

path, Z e ias — an i , and to an elastic scattering time 
T i as — (vfTh) • Hence, the regions with few impurities 
can be considered low-density metals in the dirty limit, 

aS T c7as - £ F- 

The Dirac equation allows for localized solutions that 
satisfy many possible boundary conditions. It is known 
that small circu lar defects result in localized and semi- 
localized states ( Dong et~al\ . Il998h . that is, states whose 



wavefunction decays as 1/r as a function of the distance 
from the center of the defect. A discrete version of 
these states can be realized in a nearest neighbor tight- 
binding model with unitary scatterers such as vacancies 
( Pereira et al, 2006). n the continuum, the Dirac equa- 
tion (fl9| for the wavefunction, ^(r) = (0i(r), </>2(r)), can 
be written as: 



d w <f>i(w,w*) 



0, 
0, 



(157) 



where w — x + iy is a complex number. These equa- 
tions indicate that at zero energy the components of the 
wavefunction can only be either holomorphic or anti- 
holomorphic with respect to the variable w (that is, 
4>\(w,w*) = (pi(w*) and <f>2(w,W*) = 4>2(w)). Since the 
boundary conditions require that the wavefunction van- 
ishes at infinity the only possible solutions have the form: 
oc (l/(x + iy) n ,0) or oc (0, l/(x - iy) n ). 

The wavefunctions in the discrete lattice must be real, 
and at large distances the actual solution found near 
a vacancy tends to a superposition of two solutions 
formed from wavefunctions from the two valleys with 
equal weight, in a way s imilar to the mixing at armchair 
edges ( Brev and Fertigl . l2006bl ). 

The construction of the semi-localized state around a 
vacancy in the honeycomb lattice can be extended to 
other discrete models which leads to the Dirac equation 
in the continuum limit. A particular case is the nearest 
neighbor square lattice with half flux per plaquette, or 
the nearest neighbor square lattice with two flavors per 
site. The latter has been extensively studied in relation 
to the effects of impurities on the electronic stru cture 
of d-wave superconductors ( Balatskv et all . (2006'), and 
numerical results are in agreement with the existence of 
this solution. As the state is localized on one sublattice 




Figure 27 (Color online) Sketch of a rough graphene surface. 
The full line gives the boundary beyond which the lattice can 
be considered undistorted. The number of mid-gap states 
changes depending on a difference in the number of removed 
sites for two sublattices. 



only, the solution can be generalized for the case of two 
vacancies. 



D. Localized states near edges, cracks, and voids 

Localized states can be defined at edges where the 
number of atoms in the two sublattices is not compen- 
sated. The number of them depend on details of the edge. 
The graphene edges can be strongly deformed, due to the 
bonding of other atoms to carbon atoms at the edges. 
These atoms should not induce states in the graphene ir 
band. In general, a boundary inside the graphene ma- 
terial will exist, as sketched in Fig. [27l beyond which 
the sp 2 hybridization is well defined. If this is the case, 
the number of mid-gap states near the edge is roughly 
proportional to the difference in sites between the two 
sublattices near this boundary. 

Along a zigzag edge there is one localized state per 
three lattice units. This implies that a precursor 
structure for localized states at the Dirac energy can 
be found in ribbons or c onstrictions of small lengths 
( Munoz-Rojas et aZl . l2006[ ). which modifies the electronic 
structure and transport properties. 

Localized solutions can also be found near other defects 
that contain broken bonds or vacancies. These states do 
not allow an analytical solution, although, as discussed 
above, the continuum Dirac equation is compatible with 
many boundary conditions, and it should describe well lo- 
calized states that vary slowly over distances comparable 
to the lattice spacing. The existence of these states can 
be investigated by analyzing the scaling of the spectrum 
near a defect as a functio n of the size of the system, L 
(jVozmediano et all . 120051 ) . A number of small voids and 
elongated cracks show states whose energy scales as L~ 2 , 
while the energy of extended states scales as L~ 1 . A state 
with energy scaling L~ 2 is compatible with continuum 
states for which the modulus of the wavefunction decays 
as r~ 2 as a function of the distance from the defect. 
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E. Self-doping 



The band structure calculations discussed in previ- 
ous sections show that the electronic structure of a sin- 
gle graphene plane is not strictly symmetrical in energy 
(Reich et all l2002h . The absence of electron-hole sym- 



metry shifts the energy of the states localized near impu- 
rities above or below the Fermi level, leading to a transfer 
of charge from/to the clean regions. Hence, the combina- 
tion of localized defects and the lack of perfect electron- 
hole symmetry around the Dirac points leads to the pos- 
sibility of self-doping, in addition to the usual scattering 
processes. 

Extended lattice defects, like edges, grain boundaries, 
or micro-cracks, are likely to induce a number of elec- 
tronic states proportional to their length, L/a, where a 
is of the order of the lattice constant. Hence, a distribu- 
tion of extended defects of length I at a distance equal 
to L itself gives rise to a concentration of L/a carriers 
per carbon in regions of size of the order of (L/a) 2 . The 
resulting system can be considered a metal with a low 
density of carriers, n carr i or oc a/L per unit cell, and an 
elastic mean free path Z c i as ~ L. Then, we obtain: 



CF 



?'F 



Tela 



(158) 



and, therefore, (r i as ) _1 <C cf when a/L <C 1. Hence, 
the existence of extended defects leads to the possibil- 
ity of self-doping but maintaining most of the sample in 
the clean limit. In this regime, coherent oscillations of 
transport properties are expected, although the observed 
electronic properties may correspond to a shifted Fermi 
energy with respect to the nominally neutral defect-free 
system. 

One can describe the effects that break electron-hole 
symmetry near the Dirac points in terms of a finite next- 
nearest neighbor hopping between ir orbitals, £', in (| 148|) . 
Consider, for instance, electronic structure of a ribbon of 
width L terminated by zigzag edges, which, as discussed, 
lead to surface states for t' = 0. The translational sym- 
metry along the axis of the ribbon allows us to define 
bands in terms of a wavevector parallel to this axis. On 
the other hand, the localized surface bands, extending 
from fc|| = (27r)/3 to k\\ = — (2ir)/3 acquire a dispersion 
of order t' . Hence, if the Fermi energy remains unchanged 
at the position of the Dirac points (eDirac = St'), this 
band will be filled, and the ribbon will no longer be charge 
neutral. In order to restore charge neutrality, the Fermi 
level needs to be shifted by an amount of the order of 
t'. As a consequence, some of the extended states near 
the Dirac points are filled, leading to the phenomenon of 
self-doping. The local charge is a function of distance to 
the edges, setting the Fermi energy so that the ribbon is 
globally neutral. Note that the charge transferred to the 
surface states is mostly localized near the edges of the 
system. 
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Figure 28 (Color online) Top: Self-consistent analysis of the 
displaced charge density (in units of number of electrons per 
carbon) is shown as a continuous line, and the corresponding 
electrostatic potential (in units of t) is shown as a dashed 
line, for a graphene ribbon with periodic boundary conditions 
along the zig-zag edge (with a length of L = 960a) and with a 
circumference of size W = 80y/3a. The inset shows the charge 
density near the edge. Due to the presence of the edge, there 
is a displaced charge in the bulk (bottom panel) that is shown 
as a function of width W. Notice that the displaced charge 
vanishes in th e bulk limit (W — > o o), in agreement with (1161 [) . 
Adapted from IPeres et all l2006d . 



The charge transfer is suppressed by electrostatic ef- 
fects, as large deviation s from charge neutr ality have an 
associated energy cost ( Peres et In order to 

study these charging effects we add to the free-electron 
Hamiltonian ((5]) the Coulomb energy of interaction be- 
tween electrons: 



H, 



(159) 



where rij = J2a(. a i CT a i,o-+&j -&j,a-) is the number operator 



at site Rj, and 



J eo|R-?: — R-i" | 



(160) 



is the Coulomb interaction between electrons. We expect, 
on physics grounds, that an electrostatic potential builds 
up at the edges, shifting the position of the surface states, 
and reducing the charge transferred to/from them. The 
potential at the edge induced by a constant doping 5 
per carbon atom is roughly, ~ (8e 2 /a)(W/a) (Se 2 /a is 
the Coulomb energy per carbon), and W the width of 
the ribbon (W/a is the number of atoms involved). The 
charge transfer is stopped when the potential shifts the 
localized states to the Fermi energy, that is, when t' « 
(e 2 /a)(W/a)S. The resulting self-doping is therefore 



fa 



(161) 
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Figure 29 (Color online) Gauge field induced by a simple elas- 
tic strain. Top: The hopping along the horizontal bonds is 
assumed to be changed on the right hand side of the graphene 
lattice, defining a straight boundary between the unperturbed 
and perturbed regions (green dashed line). Bottom: The 
modified hopping acts like a constant gauge field, which dis- 
places the Dirac cones in opposite directions at the K and K' 
points of the Brillouin zone. The conservation of energy and 
momentum parallel to the boundary leads to a deflection of 
electrons by the boundary. 



that vanishes when W — * oo. 

We treat Hamiltonian (|159|1 within the Hartree approx- 
imation (that is, we replace Hi by Hm.f. = J2i Vi n i 
where Vi = YljUi,j(nj), and solve the problem self- 
consistently for (rij}). Numerical results for graphene 
ribbons of length L = 80\/3a and different widths are 



shown in Fig. [28] (t'/t = 0.2 and e 2 /a = 0.5t). The 
largest width studied is ~ 0.1/^m, and the total number 
of carbon atoms in the ribbon is w 10 5 . Notice that as 
W increases, the self-doping decreases indicating that, 
for a perfect graphene plane (W — > oo), the self-doping 
effect disappears. For realistic parameters, we find that 
the amount of self-doping is 10~ 4 — 10~ 5 electrons per 
unit cell for sizes 0.1 — 1/im. 



F. Vector potential and gauge field disorder 



As discussed in Sec. IIV1 lattice distortions modify the 
Dirac equation that describes the low energy band struc- 
ture of graphene. We consider here deformations that 
change slowly on the lattice scale, so that they do not mix 
the two inequivalent valleys. As shown earlier, perturba- 
tions that hybridize the two sublattices lead to terms 
that change the Dirac Hamiltonian from vpcr ■ k into 
vper ■ k + er • A. Hence, the vector A can be thought of as 
if induced by an effective gauge field, A. In order to pre- 
serve time reversal symmetry, this gauge field must have 
opposite signs at the two Dirac cones, Ax = —Ak 1 - 

A simple example is a distortion that changes the hop- 
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Table II Estimates of the effective magnetic length, and effec- 
tive magnetic fields generated by the deformations considered 
in this section. The intrinsic curvature entry also refers to the 
contribution from topological defects. 



ping between all bonds along a given axis of the lattice. 
Let us assume that the sites at the ends of those bonds de- 
fine the unit cell, as sketched in Fig. [29l If the distortion 
is constant, its only effect is to displace the Dirac points 
away from the BZ corners. The two inequivalent points 
are displaced in opposite directions. This uniform dis- 
tortion is the equivalent of a constant gauge field, which 
does not change the electronic spectrum. The situation 
changes if one considers a boundary that separates two 
domains where the magnitude of the distortion is differ- 
ent. The shift of the Dirac points leads to a deflection 
of the electronic trajectories that cross the boundary, as 
also sketched in Fig. EH The modulation of the gauge 
field leads to an effective magnetic field, which is of op- 
posite sign for the two valleys. 

We have shown in Section IIV.BI how topological lat- 
tice defects, such as disclinations and dislocations, can be 
described by an effective gauge field. Those defects can 
only exist in graphene sheets that are intrinsically curved, 
and the gauge field only depends on topology of the lat- 
tice. Changes in the nearest neighbor hopping also lead 
to effective gauge fields. We consider next two physical 
processes that induce effective gauge fields: (t) changes 
in the hopping induced by hybridization between it and 
a bands, which arise in curved sheets, and (ii) changes in 
the hopping due to modulation in the bond length, which 
is associated with elastic strain. The strength of these 
fields depends on parameters that describe the value of 
the ir-a hybridization, and the dependence of hopping on 
the bond length. 

A comparison of the relative strengths of the gauge 
fields induced by intrinsic curvature, it — a hybridization 
(extrinsic curvature), and elastic strains, arising from a 
ripple of typical height and size is given in Table HT1 



1. Gauge field induced by curvature 



As we discussed in Sec. IIV.A} when the n orbitals are 
not parallel, the hybridization between them depends on 
their relative orientation. The angle 6i determines the 
relative orientation of neighboring orbitals at some posi- 
tion Yi in the graphene sheet. The value of 0j depends 
on the local curvature of the layer. The relative angle of 
rotation of two p z orbitals at position and Tj can be 
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written as: cos(9i — 6j) = N, ■ Nj, where Nj is the unit 
vector perpendicular to the surface, defined in (|125fl . If 
Tj = r, : + we can write: 



N r N i «l+N i -[(u ii -V)N i ] + -N i -[(u i ,-V) 2 N i ] ,(162) 



where we assume smoothly varying N(r). We use (|125p 
in terms of the height field h(r) (N(r) m z — V/i(r) — 
(Vh) 2 z/2) to rewrite fTM]) as: 



N ! .N J «l--[(uyV)V/i(r 1 



(163) 



Hence, bending of the graphene sheet leads to a modifi- 
cation of the hopping amplitude between different sites 
in the form: 



%~-y[K- -V)V/i(r 4 )] 2 , 



(164) 



where ^ is the bare hopping energy. A similar ef- 
fect leads t o changes the e l ectron ic states in a carbon 
nanotubes ( Kane and Meld . Il997j) . Using the results 
of Sec. llVl namely l|147j) . we can now see that a vec- 
tor potent ial is generated for nearest neighbor hopping 
(u = Sab) (|Eun-Ah Kim and Castro Netd . l2007j) : 



Al h) = [{dlhf - (d 2 y hf 



8 

3E ab a 2 



{d 2 x h + d 2 y h) d x hd y h (165) 



where t he coupling constant E n b depends on m icroscopic 
details l|Eun-Ah Kim and Castro Netd . [2007h . The flux 
of effective magnetic field through a ripple of lateral di- 
mension I and height h is given approximately by: 



2(,2 



E ab a 2 h 
I 3 



(166) 



where the radius of curvature is R 1 ~ hi 2 . For a ripple 
with I w 20nm, h ~ lnm, taking E a b/vp ~ 10 A" 1 , we 
find $ w 10- 3 $ - 



2. Elastic strain 

The elastic free energy for graphene can be written in 
terms of the in-plane displacement u(r) = (u x , u y ) as: 



F[u] = ±/ 



i=l,2 



(167) 



There are many types of static deformation of the hon- 
eycomb lattice which can affect the propagation of Dirac 
fermions. The simplest one is due to changes in the 
area of the unit cell either due to dilation or contrac- 
tion. Changes in the unit cell area lead to local changes 
in the density of electrons and, therefore, local changes in 
the chemical potential in the system. In this case, their 
effect is similar to the one found in (| 148[) . and we must 
have: 



dp (r) = g{u xx + u yy ) , 



(169) 



and their effect is diagonal in the sublattice index. 

The nearest neighbor hopping depends on the length of 
the carbon bond. Hence, elastic strains that modify the 
relative orientation of the atoms also lead to an effective 
gauge field, which acts on each K point separ ately, as 
first discussed in relation to car bon nanotubes I Manesl . 
120071 : IS uzuura and Andd . r2002bl ). Consider two carbon 
atoms located in two different sublattices in the same 
unit cell at R;. The change in the local bond length can 
be written as: 

Sm = — ■ [u A (Ri) - u B (R 4 + S ab )} . (170) 
a 

The local displacemen ts of the atom s in the unit cell can 
be related to u(r) by (jAndd . l2006ah : 



(S ab ■ V)u = k 1 (u A - u B ) 



(171) 



where k is a dimensionless quantity that depends on mi- 
croscopic details. Changes in the bond length lead to 
changes in the hopping amplitude: 



U 



Of 



da 



and we can write: 



where 



St {ab) ( 



Su(r) 



Q t (ab) 
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(172) 



(173) 



(174) 



Substituting (1 1701) int o HI 731) and the final result into 
(fT47| . one finds (|Andd . l2006al ): 



3 

tA.y — — j3 H Uxy ■ 



(175) 



where B is the bulk modulus, Q is the shear modulus, 
and 



it 



1 / dui duj 



is the strain tensor (xi = x and x-i = y), 



We assume that the strains induced by a ripple of di- 
mension I and height h scale as Uij ~ (h/l) 2 - Then, 



using (3/vy 



lA 1 , we find that the total flux 



(168) through a ripple is: 



__ 

al 



(176) 
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For ripples such that h <■ 
timate gives <f> ~ 10 _1 $c 



' lnm and I ~ 20nm, this es- 
in reasonable agreement with 



the estimates in ref. (jMorozov et al\ . 12006 

The strain tensor must satisfy some additional con- 
straints, as it is derived from a displacement vector field. 
These cons traints are called Saint Ven ant compatibility 
conditions 1 Landau and Lifshitd . fl959h : 



dui 



ijkl 



+ 



dxkdxi dxidx 



dug 
dxjdxk 



du 



dxidxi 



= 0. 
(177) 

An elastic deformation changes the distances in the crys- 
tal lattice and can be considered as a change in the met- 
ric: 



9ij 



(178) 



The compatibility equations l|177j) are equivalent to the 
condition that the curvature tensor derived from l|178p is 
zero. Hence, a purely elastic deformation cannot induce 
intrinsic curvature in the sheet, which only arises from 
topological defects. The eff ective fields associate d with 
elastic strains can be large 1 Morozov et all l2006t ) . lead- 



ing to significant changes in the electronic wavefunctions. 
An analysis of the resulting state, and the possible insta- 
biliti es that may occur can be found in (jGuinea et all 
l2007t ). 



3. Random gauge fields 

The preceding discussion suggests that the effective 
fields associated with lattice defects can modify signif- 
icantly the electronic properties. This is the case when 
the fields do not change appreciably on scales comparable 
to the (effective) magnetic length. The general problem 
of random gauge fields for Dirac fermions has been exten- 
sively analyzed before the current interes t in graphene 
as th e topic is also relevant for the IQH E dLudwig et~al. 
1994) and d-wave superconductivity (|Nersesvan et al 
19941 ). The one electron nature of this two dimensional 



problem makes it possible, at the Dirac energy, to map 
it onto models of interacting electrons in one dimension , 
wher e many exact results can be obtained (jCastillo et all 
Il997h . The low energy density of states, p(oj), acquires 
an anomalous exponent, p(uj) oc |w| 1_A , where A > 0. 
The density of states is enhanced near the Dirac en- 
ergy, reflecting the tendency of disorder to close gaps. 
For sufficiently large values of the random gauge field , 
a phase transition is also po ssible ( Chamon et aH , [l99ffl 
iHorovitz and Doussall . l2002h . 

Perturbation theory shows that random gauge fields 
are a marginal perturbation at the Dirac point, leading 
to logarithmic divergences. These divergences tend to 
have the opposite sign with respect to those induced by 
the Coulomb interaction (see Sec. IV.Bjl . As a result, a 
renormalization group (RG) analysis of interacting elec- 
trons in a random g auge field suggests the possibility 
of non-trivial phases ijAleiner and Efetovl . l2006t lAltlandl . 
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G. Coupling to magnetic impurities 

Magnetic impurities in graphene can be intro- 
duced chemically b y depo siti on and int erc alatio n 
(|Calandra and Mauril . l2007t lUchoa et all 120071 ). 
or self-generated by th e intr o ductio n of defects 
(Kumazaki and Hirashimal . 120061 . l2007t ). The en- 
ergy dependence of the density of states in graphene 
leads to changes in the formation of a Kondo reso- 
nance between a magnetic impurity and the graphene 
electrons. The vanishing of the density of states at 
the Dirac energy implies that a Kondo singlet in the 
ground state is not formed unless the exchange coupling 
exceeds a critical value, of the order of the electron 
bandwidth, a problem already studied in connection 
with magnetic impuriti e s in d-waye superconductors 
dCassanello and Fradkml Il996l Il997t iFritz etail . 120061 : 
IPolkovnikovL l2002t IPolkovnikov et al\ . l200lh . For weak 
exchange couplings, the magnetic impurity remains 
unscreened. An external gate changes the chemical 
potential, allowing for a t uning of the Kondo resonance 
( Sengupta and Baskaranl . l2007h . The situation changes 
significantly if the scalar potential induced by the 
magnetic impurity is taken into account. This potential 
that can be comparable to the bandwidth allows the 
formation of mid-gap states and changes the phase-shift 
assoc iated to spin scattering ( Hentschel and Guinea! . 
l2007f ). These phase-shifts have a weak logarithmic 
dependence on the chemical potential, and a Kondo 
resonance can exist, even close to the Dirac energy. 

The RKKY interaction between magnetic impurities is 
also modified in graphene. At finite fillings, the absence 
of intra-valley backscattering leads to a reduction of the 
Friede l oscillations, which decay as sin(2fcF7")/H' 
l2006bl : ICheianov and Fal'kol . 120061 : IWunsch et al. 
This effect leads to an RKKY interaction, at finite fill- 
ings, which oscillate and decay as \r\~ 3 . When interval- 
ley scattering is included, the interaction reverts to the 
usual dependence on distan ce in two dimensions, |r| -2 
(|Cheianov and FaPkol . f2006h . At half-filling extended de- 
fects lead to an RKKY interaction with an H~ 3 de - 
pendence i Dugaev et a/J . l2~006l : IVozmediano et a/1 120051 ). 




This behavior is changed wh en the impurity potential 



is loc alized on atomic scales I Brey et all . l2007t 
2007), or for highly symmetrical couplings 



2007). 
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H. Weak and strong localization 

In sufficiently clean systems, where the Fermi wave- 
length is much shorter than the mean free path, k F l ^> 1, 
electronic transport can be described in classical terms, 
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assuming that electrons follow well defined trajectories. 
At low temperatures, when electrons remain coherent 
over long distances, quantum effects lead to interference 
corrections to the classical expressions for the conduc- 
tivity, the weak localizat i on cor rection ( Bergman! . 1 19841 : 
IChakravarty and Schmidl . [l986). These corrections are 
usually due to the positive interference between two paths 
along closed loops, traversed in opposite directions. As a 
result, the probability that the electron goes back to the 
origin is enhanced, so that quantum corrections decrease 
the conductivity. These interferences are suppressed for 
paths longer than the dephasing length, 1$, determined 
by interactions between the electron and environment. 
Interference effects can also be suppressed by magnetic 
fields that break down time reversal symmetry and adds 
a random relative phase to the process discussed above. 
Hence, in most metals, the conductivity increases when 
a small magnetic field is applied (negative magnetoresis- 
tance). 

Graphene is special in this respect, due to the chiral- 
ity of its electrons. The motion along a closed path in- 
duces a change in the relative weight of the two compo- 
nents of the wavefunction, leading to a new phase, which 
contributes to the interference processes. If the electron 
traverses a path without being scattered from one val- 
ley to the other, this (Berry) phase changes the sign 
of the amplitude of one path with respect to the time- 
reversed path. As a consequence, the two paths inter- 
fere destructively, leading to a su ppression of backscatter- 
ing I Suzuura and Andol . 2002al ). Similar processes take 
place in materials with strong spin orbit coupling, as 
the s pin direction changes along the path of t h e elec - 
tron ijBergmanl . Il984t IChakravarty and Schmidl . 1 19861) . 
Hence, if scattering between valleys in graphene can be 
neglected, one expects a positive magnetoresistance, i. 
e., weak anti- localization. In general, intra- and interval- 
ley elastic scattering can be described in terms of two 
different scattering times, Ti ntra and Ti n t er , so that if 
Tintra "C Tinter one expects weak anti-localization pro- 
cesses, while if Ti n t er *C Tmtra ordinary weak localiza- 
tion will take place. Experimentally, localization ef- 
fects are always strongly suppressed close to the Dirac 
point but can be partially or, in rare cases, completely 
recovered at high carrier concent rations, depend i ng on 
a particular single-layer sample ( Morozov et all l2006t 
iTikhonenko et all l2007h . Multilayer samples exhibit an 
additional positive magnetoresistance in higher magnetic 
fields, which can be attribued to classical changes in the 
current dis tribution due to a v ertical gradient of con- 
centration ( Morozo v et al ]. l2006i ) and anti-localization ef- 
fects (|Wu et aZ.l . l2007h . 

The propagation of an electron in the absence of in- 
tervalley scattering can be affected by the effective gauge 
fields induced by lattice defects and curvature. These 
fields can s uppress the interferen c e corrections to the con- 
ducti vity I Morozov et all 120061 : iMorpurgo and Guineal . 
l2006h . In addition, the description in terms of free Dirac 
electrons is only valid near the neutrality point. The 



Fermi energy acquires a trigonal distortion away from 
the Dirac point, and backward scatteri ng within each 
valley is no longer completely suppressed ( McCann et all 
leading to a further suppression of anti-localization 
effects at high dopings. Finally, the gradient of external 
potentials i nduce a small asymme t ry be tween the two 
sublattices ( Morpurgo and Guineal . 120061 ). This effect 
will also contribute to reduce anti-localization, without 
giving rise to localization effects. 

The above analysis has to be modified for a graphene 
bilayer. Although the description of the electronic states 
requires a two component spinor, the total phase around 
a closed loop is 2 i r and backscattering is not suppressed 
( Kechedzhi et all l2007h . This result is consistent with 
experimental observations, which sho w the existence of 
weak localization effects in a bilayer ( Gorbachev et all 
l20Q7h . 

When the Fermi energy is at the Dirac point, a replica 
analysis shows that the conductivit y approaches a u ni- 
versal value of the order of e 2 /h ijFradkinl . Il986af bl). 
This res ult is valid when intervalley scattering is ne - 
glected ijOstrovskv et all 120061 . l2007t iRvu et all 120071 ) . 
Localization is induced when these term s are included 
(jAleiner and Efetovl . 120061 : lAltlaiidl. l2006h. a s also con- 
firmed by numerical calculations ()Louis et «/.l . l2007f l. In- 
teraction effects tend to suppress the effects of disorder. 
The same result, namely a conductance of the order of 
e 2 /h, is obtained for disordered graphene bilayers where a 
self-consistent calculati on leads to universal conductivity 
at the neutral ity point ( Katsnelsonl . l2007cl : lNilsson et all 
2006a, 2007a! ) . In a biased graphene bilayer, the pres- 
ence of impurities leads to the appearance of impurity 
tails in the density of states due to the creation of mid- 
gap states which are sensitive to the applied electric field 
that o pens the gap between the condu ction and valence 
bands (jNilsson and Castro Netol . l2007t ). 

One should point out that most of the calculations 
of transport properties assume self-averaging, that is, 
that one can exchange a problem with lack of transla- 
tional invariance by an effective medium system with 
damping. Obviously this procedure only works when 
the disorder is weak and the system is in the metallic 
phase. Close to the localized phase this procedure breaks 
down, the system divides itself into regions of different 
chemical potential and one has to think about transport 
i n real space, usually d e scribed in t erms of percolation 
()Cheianov et all l2007bt IShklovskiil . I2007L Single elec- 
tron transistor (SET) measurements of graphene show 
that t his seems to be the situation in graphene at half- 
filling (jMartin et all l2007tl . 

Finally, we should point out that graphene stacks 
suffer from another source of disorder, namely, c-axis 
disorder that is either due to impurities between lay- 
ers or rotation of graphene planes relative to each 
other. In either case the in-plane and out-of-plane 
transport is directly affected. This kind of disorder 
has been observed experiment a lly by d ifferent techniques 
()Bar et all l2007t lHass et all l2007bh . In the case of 



the bilayer, the rotation of planes changes substantially 
the spectrum restoring the D irac fermion description 
( Lopes dos Santos et all l2007t ). The transport proper- 
ties in the out of plane direction are determined by the in- 
terlayer current operator, j„, n+ i = itY!,( c A,n,s c An+l,s - 

c \ n+i s c A.n,s), where n is a layer index, and A is a 
generic index that defines the sites coupled by the in- 
terlayer hopping t. If we only consider hopping between 
nearest neighbor sites in consecutive layers, these sites 
belong to one of the two sublattices in each layer. 

In a multilayer with Bernal stacking, these connected 
sites are the ones where the density of states vanishes at 
zero energy, as discussed above. Hence, even in a clean 
system, the number of conducting channels in the direc- 
tion perpendicular to the lay ers vanishes at zero energy 
(jNilsson et all l2006aL l2007aT) . This situation is reminis- 
cent of the in plane transport properties of a single layer 
graphene. Similar to the latter self-consistent 
Born approximation for a small concentration of impu- 
rities leads to a finite conductivity, which becomes inde- 
pendent of the number of impurities. 



I. Transport near the Dirac point 

In clean graphene, the number of channels available 
for electron transport decreases as the chemical poten- 
tial approaches the Dirac energy. As a result, the con- 
ductance through a clean graphene ribbon is, at most, 
4e 2 //i, where the factor of 4 stands for the spin and val- 
ley degeneracy. In addition, only one out of every three 
possible clean graphene ribbons have a conduction chan- 
nel at the Dirac energy. The other two thirds are semi- 
conducting, with a gap of the order of vf/W, where W 
is the width. This result is a consequence of the addi- 
tional periodicity introduced by the wavefunctions at the 
K and K' points of the Brillouin Zone, irrespective of 
the boundary conditions. 

A wide graphene ribbon allows for many channels, 
which can be approximately classified by the momen- 
tum perpendicular to the axis of the ribbon, k y . At the 
Dirac energy, transport through these channels is inhib- 
ited by the existence of a gap, Ak y = vpk y . Transport 
through these channels is suppressed by a factor of the 
order of e~ kyL , where L is the length of the ribbon. The 
number of transverse channels increases as W/a, where 
W is the width of the ribbon and a is a length of the 
order of the lattice spacing. The allowed values of k y 
are oc n y /W, where n y is an integer. Hence, for a rib- 
bon such that W 3> L, there are many channels which 
satisfy k y L <C 1. Transport through these channels is 
not strongly inhibited, and their contribution dominates 
when the Fermi energy lies near the Dirac point. The 
conductanc e arising from these channels is given appro x- 
imately by ijKatsnelsonl . l2006bl iTworzvdlo et all 12006? ) : 



The transmission at normal incidence, k y = 0, is one, 
in agreement with the absence of backscattering in 
graphene, for any barrier that does not induce intervalley 
scattering ( Katsnelson et «ill2006f ). The transmission of 
a given channel scales as T(k y ) = 1/ cosh 2 (fcj,L/2). 

Eq. (|179p shows that the contribution from all trans- 
verse channels lead to a conductance which scales, simi- 
lar to a function of the length and width of the system, 
as the conductivity of a diffusive metal. Moreover, the 
value of the effective conductivity is of the order of e 2 /h. 
It can also be shown that the shot noise depends on cur- 
rent in the same way as in a diffusive metal. A detailed 
analysis of possible boundary conditions at the contacts 
and their influence on evanescent waves can be fou nd in 
(|Robinson and Schomerul 120071 : ISchomerusl l2007h . The 
calculations leading to eq. (|179[) can be extended to a 
graphene bilayer. The conductance is, again, a sum- 
mation of terms arising from evanescent waves between 
the two contacts, and it has the dependence on sam- 
ple dimensions of a 2D conduc tivity of the order of e 2 /h 
( Snvman and Beenakkerll2007l ). although there is a pref- 
actor twice bigger than the one in single layer graphene. 

The calculation of the conductance of clean graphene 
in terms of transmission coefficients, using the Landauer 
method leads to an effective conductivity which is equal 
to the value obtained for bulk grap hene using diagram - 
matic methods, the Kubo formula 1 Peres e^oll . l2006dh . 
in the limit of zero impurity concentration and zero dop- 
ing. Moreover, this correspondence remains valid for the 
case of a bilayer without and with trigonal warp ing effects 
()Cserti et all l2007al : iKoshino and Andcl l200rl ) . 

Disorder at the Dirac energy changes the conduc- 
tance of graphene ribbons in two opposite directions 
( Louis etal\ . 120071 ): i) a sufficiently strong disorder, 
with short range (intervalley) contributions, lead to a 
localized regime, where the conductance depends ex- 
ponentially on the ribbon length, and ii) at the Dirac 
energy, disorder allows mid-gap states that can en- 
hance the conductance mediated by evanescent waves 
discussed above. A fluctuating electrostatic potential 
also reduces the effective gap for the transverse chan- 
nels, enhancing further the conductance. The resonant 
tunneling regime mediated by mi d-gap s t ate w as sug- 
gested by analytical calculations I Titovl . l2007h . The 
enhancement of the conductance by potential fluctu- 
ations can also be studied semi-analytically. In the 
absence of intervalley scattering, it leads to an ef- 
fective conduct i vity which grows with ribbon length 
(ISan-Jose et all 120071). In fact, analytical and numerical 
studies (Bardarson et all 2007 : Lewenkqpf et all . 120071 : 



iNomura et all 2007 ; San- Jose et al. , 2007 ) show that the 
conductivity obeys a universal scaling with the lattice 
size L: 



*(L) = ^-(Aln(L/Z) + B) 



(180) 



h 2tt 



— fay L 



h L 



(179) 



where £ is a length scale associated with range of inter- 
actions and A and B are numbers of the order of unit 
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(A m 0.17 and B ps 0.23 for a graphene lattice in the 
shape of a square of size L (|Lewenkopf et aD . l2007h ) . No- 
tice, therefore, that the conductivity is always of the or- 
der of e 2 jh and has a weak dependence on size. 



J. Boltzmann Equation description of DC transport in doped 
graphene 



It was shown experimentally that the DC conductiv- 
i ty of graphene depend s linearly on the gate potential 
(|Novoselov et all l2005aL 120041 l2005bh . except very close 
to the neutrality point (see Fig l30|) . Since the gate po- 
tential depends linearly on the electronic density, n, one 
has a conductivity a oc n. As shown by Shon and Ando 
( Shon and Andol . ll998l ) if the scatterers are short range 



The Boltzmann equation has the form 1 Zimanl . Il972h 



the DC conductivity should be independent of the elec- 
tronic density, a t odds with the experimental result. It 
has b e en sh own 1 Andol . [2006bl : iNomura and MacDonaldl . 
20061 . l2007h that by considering a scattering mech- 
anism based on screened charged impurities it is 
possible to obtain from a Boltzmann equation ap- 
proach a conductivity varying linearly with the den- 
sity, in agreement with the e x perimental resul t 
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Figure 30 (Color online) An example of changes in conduc- 
tivity a of graphene with varying gate voltage, V g , and, there- 
fore, carrier concentration n. Here a is proportional to n as 
discussed in the text. Note that samples with higher mo- 
bility (> 1 m 2 /Vs) normally show a sublinear dependence, 
presumably indicating the presence of different types of scat- 
terers. Inset: scanning-electron micrograph of one of experi- 
mental devices (in false colors matching those seen in visible 
optics. The scale of the micrograph is given by the width 
of the Hall bar, which is one micrometer. Adapted from 
ijNovoselov et oil 120053 ). 



-v k -V r f(e k )-e(E+v k xH)-V k f(e k ) = - ^ 

^ scatt. 

(181) 

The solution of the Boltzmann equation in its general 
form is difficult and one needs therefore to rely upon 
some approximation. The first step in the usual approx- 
imation scheme is to write the distribution as f(e k ) = 
f°i e k) + <?(efc) where f°(e k ) is the steady state distri- 
bution function and g(e k ) is assumed to be small. In- 
serting this ansatz in (|18ip and keeping only terms that 
are linear in the exter nal fields one obtains the linearized 
Boltzmann equation ( Zimanl . Il972f ) which reads 



df°(e k ) 
de k 

dfk 
dt 



v k 



V r T + e( E -\v r Q 



v k ■ V r g k + e(v k x H) ■ V k gk ■ 



(182) 



The second approximation has to do with the form of the 
scattering term. The simplest approach is to introduce a 
relaxation time approximation: 



dt 



scatt. ^ 



(183) 



where Tk is the relaxation time, assumed to be momen- 
tum dependent. This momentum dependence is deter- 
mined phenomenologically in such way that the depen- 
dence of the conductivity upon the electronic density 
agrees with experimental data. The Boltzmann equation 
is certainly not valid at the Dirac point, but since many 
experiments are performed at finite carrier density, con- 
trolled by an external gate voltage, we expect the Boltz- 
mann equation to g ive reliable results if an appropriate 
form for Tk is used ( Adam et all . 120071 ) . 

Let us compute the Boltzmann relaxation time, Tk, 
for two different scattering potentials: (i) a Dirac delta 
function potential; (ii) a unscreened Coulomb potential. 
The relaxation time T k is defined as: 

— =th [de /^S(fc,fc')(l-cos0), (184) 

Tk J J \^) 

where m is impurity concentration per unit of area, and 
the transition rate S(k, k') is given, in the Born approx- 
imation, by: 

S(k,k') = 2TT\H k , tk \ 2 — 5(k'-k), (185) 
vf 

where the Vffc is the dispersion of Dirac fermions in 
graphene and Hk\k is defined as 



(186) 



H k ',k= / drip%,(r)Us(r)ipit(r) , 



with [/g(r) the scattering potential and tp k (r) is the elec- 
tronic spinor wavefunction of a clean graphene sheet. If 



the potential is short range. i Shon and Andol . fT998f ) of the 
form Us = vo6(r), the Boltzmann relaxation time turns 
out to be 



Tfc 



Av F 1 
mv 2 k 



(187) 



On the other hand, if the potential is the Coulomb po- 
tential, given by Us(r) = eQ / (Aireoer) for charged impu- 
rities of charge Q, the relaxation time is given by 



Th 



(188) 



where Uq 



riiQ 2 e 2 / (lQe 2 e 2 ) . As we argue below, the 
phenomenology of Dirac fermions implies that the scat- 
tering in graphene must be of the form (|188|l . 

Within the relaxation time approximation the solution 
of the linearized Boltzmann equation when an electric 
field is applied to the sample is 



9k 



df°(e k ) 



er k v k ■ E , 



(189) 



and the electric current reads (spin and valley indexes 
included) 



J = ev k9k ■ 



(190) 
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The problem of magnetotransport in the presence of 
Coulomb impurities, as discussed in the previous sec- 
tion is still an open research problem. A similar ap- 
proach was considered by Abrikosov in the quantum mag- 
netoresistance stu dy of non-stoichiometric chalcogenides 
(jAbrikosovl Il998t ). In the case of graphene, the effec- 
tive Hamiltonian describing Dirac fermions in a magnetic 
field (including disorder) can be written as: H = H + Hi 
where H n is given by OP) and Hi is the impurity potential 
reading 1 Peres et all . l2006c 




(192) 



The formulation of the problem in second quantization 
requires the solution of Hq, which was done in Section 
III. II The field operators, close to the K point, are defined 
as (the spin index is omitted for simplicity): 



*(r) = 




E 



. ik.r 



Ck-l 



n (y - kl 2 B ) 



2L \ (j) n+ i(y~ kl 2 B ) 



,(193) 



Since at low temperatures the following relation 
—df°(e k )/de k — > S(fi — vpk) holds, one can easily see 
that assuming l|188j) where k is measured relatively to 
the Dirac point, the electronic conductivity turns out to 
be 



2 £V 

h Uq 



2 9 



(191) 



where uq is the strength of the scattering potential 
(with dimensions of energy). The electronic conductiv- 
ity depends linearly on the electron density, in agree- 
ment with the experimental data. We stress that the 
Coulomb potential is one possible mechanism of produc- 
ing a scattering rate of the form (|188|l but we do not ex- 
clude that other mechanisms may exist (see, for instance, 
()Katsnelson and Geiml . l2008h ). 



K. Magnetotransport and universal conductivity 

The description of the magnetotransport properties of 
electrons in a disordered honeycomb lattice is complex 
because of the inte rference effects associated with the 
Hofstadter problem i Gumbs and Feketel . ll997l ). We shall 
simplify our problem by describing electrons in the hon- 
eycomb lattice as Dirac fermions in the continuum ap- 
proximation, introduced in Sec. III.BI Furthermore, we 
will only focus on the problem of short range scatter- 
ing in the unitary limit sinc e in this regime many an- 
alytical results are obtained ( Kumazaki and Hirashimal . 



where Ck } n,a destroys an electron in band a = ±1, with 
energy level n and guiding center kl B ; c k -i destroys an 
electron in the zero Landau level; the cyclotron frequency 
is given by (|96|l . The sum over n = 0, 1, 2, . . . , is cut off at 
no given by E(l, no) = W, where W is of the order of the 
electronic bandwidth. In this representation Hq becomes 
diagonal, leading to Green's functions of the form (in 
Matsubara representation): 



G*o(fc, n, a; iu>) 



1 



iuj — E(a, n) ' 



where 



E(a, n) = auj c 



(194) 



(195) 



are the Landau levels for this problem (a = ±1 labels the 
two bands). Notice that Go(k, n, a; iuj) is effectively k- 
independent, and E(a, —1) = is the zero energy Landau 
level. When expressed in the Landau basis, the scatter- 
ing Hamiltonian l|192p connects Landau levels of negative 
and positive energy. 



1. The full self-consistent Born approximation (FSBA) 

In order to describe the effect of impurity scattering 
on the magnetoresistance of graphene, the Green's func- 
tion for Landau levels in the presence of disorder needs 
to be computed. In the context of the 2D electron gas, 
an equivalent study was performed by Ohta and Ando, 
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Figure 31 (Color online) Top: Electronic density of states 
(DOS), p(uj), as a function of lu/lo c (lu c = 0.14 eV) in a mag- 
netic field B = 12 T for different impurity concentrations n;. 
Bottom: p(ui), as a function of w/w c (lj c = 0.1 eV is the cy- 
clotron frequency) in a magnetic field B = 6 T. The solid line 
shows the DOS in the absence of disorder. The position of 
the Landau levels in the absence of disorder are shown as ver- 
tical lines. The two arrows in the top panel show the position 
of the renormalized Landau levels f see Fig|32l given b y the 
solution of Eq. (202J. Adapted from lPeres et aZ.l l2006d . 



(lAndol. Il974al lbllcl. Il975t lAndo and Uemural . 1 1974 lOhtal . 
Il968l . ll9Tlh using the averagin g procedure over impurity 
positions of Duke l Dukel . fT968t ). Below the averaging pro- 
cedure over impurity positions is performed in the stan- 
dard way, namely, having determined the Green's func- 
tion for a given impurity configuration (ri, . . . rjv 4 ), the 
position averaged Green's function is determined from: 

(G(p,n,a;iu);n, ...r Ni ))= G(p,n,a;iu;) 



L -2N t 



n/^ 

3 = 1 J 



G(p,n,a;itu]r 1 , . . .rjy 4 ) • (196) 



In the presence of Landau levels the average over im- 
purity positions involves the wavefunctions of the one- 
dimensional harmonic oscillator. After a lengthy algebra, 
the Green's function in the presence of vacancies, in the 
FSBA, can be written as: 

G(p,n,a;uj + + ) = [u - E(n, a) - Ei(w)]- 1 (197) 
G(p, -l;u + + ) = [u - M")}' 1 , (198) 



where 



Ei(w) = -mlZ^)}- 1 , (199) 

E 2 (w) = -rii[g c G{p,-l;oj + 0+)/2 + Z (tu)}- 1 (200) 

Z(u) = g c G(p,-l;u + + )/2 

+ g c ^2G(p,n,a;cj + + )/2, (201) 



and g c = A c /{2itl 2 B ) is the degeneracy of a Landau level 
per unit cell. One should notice that the Green's func- 
tions do not depend upon p explicitly. The self-consistent 
solution of Eqs. lfT97|l . ifTM]) . (p9)l . PUO]) and (f20T|) gives 
the density of states, the electron self-energy, and the 
change of Landau level energy position due to disorder. 

The effect of disorder in the density of states of Dirac 
fermions in a magnetic field is shown in Fig. [3TJ For ref- 
erence we note that £7(1, 1) = 0.14 eV, for B = 14 T, and 
£7(1, 1) = 0.1 eV, for B = 6 T. From Fig. EH we see that, 
for a given n, , the effect of broadening due to impurities 
is less effective as the magnetic field increases. It is also 
clear that the position of Landau levels is renormalized 
relatively to the non-disordered case. The renormaliza- 
tion of the Landau level position can be determined from 
poles of (fT97|) and (155)1 : 



lo - E(a,ri) - ReE(a;) = 0. 



(202) 



Of course, due to the importance of scattering at low 
energies, the solution to Eq. (|202[) does not represent 
exact eigenstates of system since the imaginary part of 
the self-energy is non- vanishing. However, these energies 
do determine the form of the density of states, as we 
discuss below. 

In Fig. [32) the graphical solution to Eq. (|202|) is given 
for two different energies (£7(— 1, n), with n = 1, 2), its is 
clear that the renormalization is important for the first 
Landau level. This result is due to the increase in scat- 
tering at low energies, which is present already in the 
case of zero magnetic field. The values of w satisfying 
Eq. I|202p show up in the density of states as the en- 
ergy values where the oscillations due to the Landau level 
quantization have a maximum. In Fig. (3T) the position 
of the renormalized Landau levels is shown in the upper 
panel (marked by two arrows) , corresponding to the bare 
energies £7(— 1, n), with n = 1,2. The importance of this 
renormalization decreases with the reduction of the num- 
ber of impurities. This is clear in Fig. l3T1 where a visible 
shift toward low energies is evident when m has a small 
10% change, from = 10~ 3 to rii — 9 x 10~ 4 . 

The study of the magnetoresistance properties of the 
system requires the calculation of the conductivity ten- 
sor. We compute the current-current correlation func- 
tion and from it the conductivity tensor is derived. The 
details of the calculations are presented in 1 Peres et all . 

If we however neglect the real part of the self- 
energy, assume for ImEi(w) = T (i = 1,2) a constant 
value, and consi der that £7(1, 1 ) ^> T, these results re- 
duce to those of i Gorbar et al\ . \2002\) . 

It is instructive to consider first the case u, T — > 0, 
which leads to (a xx (0, 0) = a ): 

"ImE 1 (0)/ImE 2 (0) - 1 



a 7r 



1 + (ImEi(0)/w c ) 2 
no + 1 



720 + 1 + (ImE^O)/^) 5 



(203) 



where we include a factor 2 due to the valley degeneracy. 
In the absence of a magnetic field (lo c — > 0) the above 
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Figure 32 (Color online) Imaginary (right) and real (left) 
parts of Ei(lj) (top) and T,2(ui) (bottom), in units of u> , as 
a function of cj/oj c . The right panels also show the intercept 
of uj — E{a, n) with ReS( ctj) as required by Eq. I|202p . Adapted 
from IPeres et all l2006cj . 



expression reduces to: 

£.1 

h 7T 



CO 



1 



[ImSi(O)] 5 



(vFA^ + pmE^O)] 2 



(204) 



where we have introduced the energy cut-off, vpA. Either 
when ImSi(0) ~ ImX^O) and tu c ImEi(O) (or uq ^ 
lm£i(0)/w c , lo c = E(0, 1) = V2v F /l 2 B ), or when Av F > 
ImEi(O), in the absence of an applied field, Eqs. i|203p 
and (|204[) reduce to: 



If! 

7T h ' 



(205) 
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I98|) . This result was obtained previously 
by Ando and collaborators using the second order 
self-consistent Born a pproximation l Ando et~al\ . 120021 : 
IShon and Andol . Il998h . 

Because the DC magnetotransport properties of 
graphene are normally measured with the possibility 
of tuning its electron ic density by a gate potential 
( Novoselov et aLl . l2004f l . it is important to compute the 
conductivity kernel, since this has direct experimental 
relevance. In the the case w^Owe write the conductiv- 
ity <J xx {0,T) as: 



e 2 f°° 
0^(0, T) = — / de 
J_ 00 



df(e) 
de 



(206) 



where the cond uctivity kernel KbU ) is given in the Ap- 
pendix of Ref. ( Peres et all l2006d ). The magnetic field 
dependence of kernel K B (e) is shown in Fig. [33l Ob- 
serve that the effect of disorder is the creation of a re- 
gion where Ks{e) remains constant before it starts to 
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Figure 33 (Color online) Conductivity kernel, K(lo) (in units 
of e 2 /(irh)), as a function of energy uj for different magnetic 
fields and for rii = 10~ 3 . The horizontal lines mark the uni- 
versal limit of the conductivity per cone, ao = 2e 2 /(nk). The 
vertical lines show the position of the Landau levels in the 
absence of disorder. Adapted from lPeres et a/.Ll2006d . 



increase in energy with superimposed oscillations com- 
ing from the Landau levels. The same effect, but with 
the absence of the oscillations, was identified at the level 
of the self-consistent density of states plotted in Fig. [31j 
Together with a xx (0,T), the Hall conductivity ^^(OjT) 
allows the calculation of the resistivity tensor (|109[) . 

Let us now focus on the optical cond uctivity, a xx (u>) 
(|Gusvnin et aZ.l . l2007l : IPeres et aLl . l2006ch . This quantity 
can be probed by reflectivity experiment s in the subter - 
ahertz to mid-infrared frequency range (jBliokhl . 120051) . 
This quantity is depicted in Fig. [34] for different mag- 
netic fields. It is clear that the first peak is controlled 
by the E(l, 1) — E(l, —1), and we have checked that it 
does not obey any particular scaling form as a function 
of uj/B. On the other hand, as the effect of scattering 
becomes less important the high energy conductivity os- 
cillations start obeying the scaling lo/^/B, as we show in 
the lower right panel of Fig. fJH 



V. MANY-BODY EFFECTS 



A. Electron-phonon interactions 



In Sec. lIV.FTTI and Sec. lIV.FT2l we discussed how static 
deformations of the graphene sheet due to bending and 
strain couple to the Dirac fermions via vector potentials. 
Just as bending has to do with the flexural modes of the 
graphene sheet (as discussed in Sec. lIIII). strain fields are 
relate d to optical and acoustic modes (|Wirtz and Rubiol . 
l2004h . Given the local displacements of the atoms in each 
sublattice, and ub, the electron-phonon coupling has 
essentially the form discussed previously for static fields. 

The coupling to acoustic modes is the most straight- 
forward one, since it already appears in the elastic the- 
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Figure 34 (Color online) Frequency dependent conductiv- 
ity per cone, o-(lo) (in units of e 2 /(nh)) at T = 10 K and 
m = 10~ 3 , as a function of the energy u> (in units of uj c ) 
for different values of the magnetic field B. The vertical ar- 
rows in the upper left panel, labeled a, b, and c, show the 
positions of the transitions between different Landau levels: 
E (l, 1)_B(-1,0), E(2, l)-£(-l,0), and £(1, 1)-E(1, -1), 
respectively. The horizontal continuous lines show the value 
of the universal conductivity. The lower right panel shows the 
conductivity for different values of magnetic field as a function 
of uj/y/~B. Adapted from lPeres et ail l2006d . 



ory. If u ac is the acoustic phonon displacement, then 
the relation between this displacement and the atom dis- 
placement is given by equation l|17ip . and its coupling 
to electrons is given by the vector potential (|175|1 in the 
Dirac equation (|150p . 

For optical modes the situation is sli ghtly different 
since the optical mode displacement is I Andol . l2006al 
l2007hh : 



U 



op 



-y=(u A -u B ), 



(207) 



that is, the bond length deformation vector. To calculate 
the coupling to the electrons we can proceed as previously 
and compute the change in the nearest neighbor hopping 
energy due to the lattice distortion through l|172p . f|173|l . 
(|170p . and (|207p . Once again the electron-phonon inter- 
action becomes a problem of the coupling of the electrons 
with a vector potential as in (|150p where the components 
of the vector potential are: 



/l( P) 
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2 a 2 ' 



(208) 



where j3 = dt/d\n(a) was defined in (| 1 T4[) . Notice that 
we can write: .A op = — \[?>j2{(5 f a 2 ) a x u op . A similar 
expression is valid close to the K' point with A replaced 
by -A. 



Optical phonons are particularly important in 
graphene research because of Raman spectroscopy. The 
latter has played a particul arly important ro le in the 
study of carbon nanotubes i Saito et all . 1 19981 ) because 
of the ID character of these systems, namely, the pres- 
ence of van Hove singularities in the ID spectrum lead 
to colossal enhancements of the Raman signal that can 
be easily detected, even for a single isolated carbon nan- 
otube. In graphene the situation is rather different since 
its 2D character leads to a much smoother density of 
states (except for the van Hove singularity at high en- 
ergies of the order of the hopping energy t sw 2.8 eV). 
Nevertheless, graphene is an open surface and hence is 
readily accessible by Raman spectroscopy. In fact, it has 
played a very important role becaus e it allows th e iden 
tification of the number of planes ( Ferrari et all. 
Graf et all 120071 : iGupta et ail . 120061 : IPisana et all 



Yan et al 



2006; 



2007 



2007f l. and the study of the optical phonon 
modes in graphene, particularly the ones in the center of 
the BZ with momentum q w 0. Similar studies have b een 



performed in graphite ribbons i Cancado et all . \2004) 
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Figure 35 (Color online) Top: The continuous line is the rel- 
ative phonon frequency shift as a function of £t/wo, an d the 
dashed line is the damping of the phonon due to electron-hole 
pair creation; Bottom: (a) Electron-hole process that leads 
to phonon softening (wo > 2fi), and (b) electron-hole process 
that leads to phonon hardening (luq < 2fi). 



Let us consider the effect of the Dirac fermions on 
the optical modes. If one treats the vector potential, 
electron-phonon coupling, l|150p and (|208[) up to second 
order perturbation theory, its main effect is the polar- 
ization of the electron system by creating electron-hole 
pairs. In the QED language, the creation electron-hole 
pairs is calle d pair (electron/anti-electron) production 
( Castro Netol . l2007t l . Pair production is equivalent to 
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a renormalization of the phonon propagator by a self- 
energy that is proportional to the polarization function 
of the Dirac fermions. 

The renormalized phonon frequency, f^o(q), is given 
bv llAndol. l2006al l2007bl: ICastro Neto and Guineal . 120071 : 



lLazzeri and Mauril . 1200a ISaha et all l2007h : 



2/J 2 

fi (q) « w — 5 — x(q,w ) 



(209) 



where wo is the bare phonon frequency, and the electron- 
phonon polarization function is given by: 



X (q,c)= E /(^ 
s. s '=±r v 



cPk /[£ s (k + q)]-/[£ s ,(k)] 



(210) 
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where E s (q) is the Dirac fermion dispersion (s = +1 for 
the upper band, and s = — 1 for the lower band), and 
f[E] is the Fermi-Dirac distribution function. For Ra- 
man spectroscopy, the response of interest is at q = 
where clearly only the interband processes such that 
ss' = — 1 (that is, processes between the lower and up- 
per cones) contribute. The electron-phonon polariza- 
tion function can be easily calculated using the linearized 
Dirac fermion dispersion (7|) and the low energy density 
of states lfl~5l) : 
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where we have introduced the cut-off momentum A (sw 
1/a) so that the integral converges in the ultraviolet. At 
zero temperature, T = 0, we have f[E) = 6(fi — E) and 
we assume electron doping, fi > 0, so that /[— E] = 1 (for 
the case of hole doping, [M < 0, is obtained by electron- 
hole symmetry). The integration in l|21ip gives: 



x(o,w ) 
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where the cut-off dependent term is a contribution com- 
ing from the occupied states in the lower n band and 
hence is independent of the value of the chemical poten- 
tial. This contribution can be fully incorporated into the 
bare value of ujo in (|209[i . Hence the relative shift in the 
phonon frequency can be written as: 



Sujo_ ^ 
where 
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+ m9 (cj /2 - Li) (213) 
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is the dimensionless electron-phonon coupling. Notice 
that (|213p has a real and imaginary part. The real part 



represents the actual shift in frequency, while the imagi- 
nary part gives the damping of the phonon mode due to 
pair production (see Fig. [35]). There is a clear change 
in behavior depending whether [i is larger or smaller 
than u>q/2. For \i < ujq/2 there is a decrease in the 
phonon frequency implying that the lattice is softening, 
while for fx > ujq/2 the lattice hardens. The interpre- 
tation for this effect is also given in Fig. [35l On the 
one hand, if the frequency of the phonon is larger than 
twice the chemical potential, real electron-hole pairs are 
produced, leading to stronger screening of the electron- 
ion interaction and hence, to a softer phonon mode. At 
the same time the phonons become damped and decay. 
On the other hand, if the frequency of the phonon is 
smaller than the twice the chemical potential, the pro- 
duction of electron-hole pairs is halted by the Pauli prin- 
ciple and only virtual excitations can be generated lead- 
ing to polarization and lattice hardening. In this case, 
there is no damping and the phonon is long lived. This 
amazing result has been observed experimentally by Ra - 
man spectroscopy ( Pisana et al\ . 120071 : lYan et q/J . l2007t l. 
Electron-phonon coupling has also been investigat ed the- 
oretically in the case of a finite magnetic field 1 Andd . 
l2007at iGoerbig et ali^mk . In this case, resonant cou- 
pling occurs due to the large degeneracy of the Landau 
levels and different Raman transitions are expected as 
compared with the zero-field case. The coupling of elec- 
trons to flexural mod es on a free standing graphene s heet 
was discussed in ref. i Mariani and von Oppenl . [2007l l. 



B. Electron-electron interactions 

Of all disciplines of condensed matter physics, the 
study of electron-electron interactions is probably one 
of the most complex since it involves the understand- 
ing of the behavior of a macroscopic number of variables. 
Hence, the problem of interacting systems is a field in 
constant motion and we shall not try to give here a com- 
prehensive survey of the problem for graphene. Instead, 
we will focus on a small number of topics that are of 
current discussion in the literature. 

Since graphene is a truly 2D system, it is informative 
to compare it with the more standard 2DEG that has 
been studied extensively in the last 25 years since the 
development of heterostructures and t he discovery of the 
quantum Hall effect (for a review, see dStonel . ll992h y At 
the simplest level, metallic systems have two main kind 
of excitations: electron-hole pairs and collective modes 
such as plasmons. 

Electron-hole pairs are incoherent excitations of the 
Fermi sea and a direct result of Pauli 's exclusion princi- 
ple: an electron inside the Fermi sea with momentum k is 
excited outside the Fermi sea to a new state with momen- 
tum k + q, leaving a hole behind. The energy associated 
with such an excitation is simply: uj = e k+q — e k and for 
states close to the Fermi surface (k f» kp) their energy 
scales linearly with the excitation momentum, u> q w vpq. 
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Figure 36 (Color online) Electron-hole continuum and collec- 
tive modes of: (a) a 2DEG; (b) undoped graphene; (c) doped 
graphene. 



In a system with non-relativistic dispersion such as nor- 
mal metals and semiconductors, the electron-hole contin- 
uum is made out of intra-band transitions only and exists 
even at zero energy since it is always possible to produce 
electron-hole pairs with arbitrarily low energy close to 
the Fermi surface, as shown in Fig. 1361 a). Besides that, 
the 2DEG can also sustain collective excitations such as 
plasmons that have dispersion: w p i asmon (q) cx ^/q, and 
exist outside the electron-hole c ontinuum at sufficiently 
long wavelengths I Shungl . [l986af h 

In systems with relativistic-like dispersion, such as 
graphene, these excitations change considerably, espe- 
cially when the Fermi energy is at the Dirac point. In this 
case the Fermi surface shrinks to a point and hence intra- 
band excitations disappear and only interband transi- 
tions between the lower and upper cones can exist (see 
Fig l367 b)). Therefore, neutral graphene has no electron- 
hole excitations at low energy, instead each electron-hole 
pair costs energy and hence the electron-hole occupies 
the upper part of the energy versus momentum diagram. 
In this case, plasmons are suppressed and no coherent 
collective excitations can exist. If the chemical poten- 
tial is moved away from the Dirac point then intra-band 
excitations are restored and the electron-hole continuum 
of graphene shares features of the 2DEG and undoped 
graphene. The full electron-hole continuum of doped 
graphene is shown in Fig. 1367c). and in this case plasmon 
modes are allowed. As the chemical potential is raised 
away from the Dirac point, graphene resembles more and 
more the 2DEG. 

These features in the elementary excitations of 



graphene reflect its screening properties as well. In 
fact, the polarization and dielectric functions of undoped 
graphene are rather different from the ones of the 2DEG 
(Lindhard function). In the random phase approxima- 
tion (RPA), the polarization function can be calculate d 
analytically l|Gonzalez et aZ.l . ll993d . ll994l ; IShungi . ll986al ): 



IL(q,w) 



(215) 



and hence, for u> > v F q the polarization function is imag- 
inary indicating the damping of electron-hole pairs. No- 
tice that the static polarization function (u = 0) van- 
ishes linearly with q, indicating the lack of screening in 
the system. This polarization function has been also 
calc ulated in the presence of a finite chemical poten- 
tial (lAndol. l2006bl: [Hwang and Das Sarmal . 120071 : IShungl . 
ll986aU bl: IWunsch et aZ.I . l2006h . 

Undoped, clean graphene is a semimetal, with a van- 
ishing density of states at the Fermi level. As a result the 
linear Fermi Thomas screening length diverges, and the 
long range Coulomb interaction is not screened. At finite 
electron density n, the Thomas-Fermi screening length 
reads : 
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where 
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is the dimensionless coupling constant in the prob- 
lem (the analogue of (|143fl in the Coulomb impu- 
rity problem). Going beyond the linear Thomas- 
Fermi regime, it has been shown that the Coulomb 
law is modified llForier et all [2007bt iKatsnelsonl l2006at 
IZhang and Foglerl . |200tT )~ 

The Dirac Hamiltonian in the presence of interactions 
can be written as: 

H = -iv F J d 2 r&(r)cr- V*(r) 

-p(r)p(r') , (218) 



where 
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p(r) = #t(r)*(r), 



(219) 



is the electronic density. Observe that Coulomb in- 
teraction, unlike in QED, is assumed to be instanta- 
neous since vf/c r* 1/300 and hence retardation ef- 
fects are very small. Moreover, the photons propagate 
in 3D space whereas the electrons are confined to the 
2D graphene sheet. Hence, the Coulomb interaction 
breaks the Lorentz invariance of the problem and makes 
the m any-body situation rath er different from the one in 
QED i Bavm and Chinl . ll976h . Furthermore, the problem 
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Figure 37 Hartree-Fock self-energy diagram which leads to a 
logarithmic renormalization of the Fermi velocity. 



depends on two parameters, v-p and e 2 /eo- Under a di- 
mensional scaling, r — > Ar, t — > Xt, "J — > A -1 "J, both pa- 
rameters remain invariant. In RG language, the Coulomb 
interaction is a marginal variable, whose strength relative 
to the kinetic energy does not change upon a change in 
scale. If the units are chosen in such a way that vf is 
dimensionless, the value of e 2 /eo will also be rendered 
dimensionless. This is the case in theories considered 
renormalizable in quantum field theory. 

The Fermi velocity in graphene is comparable to that 
in half-filled metals. In solids with lattice constant a, the 
total kinetic energy per site, l/(ma 2 ), where m is the 
bare mass of the electron, is of the same order of mag- 
nitude as the electrostatic energy, e 2 /(eoa). The Fermi 
velocity for fillings of the order of unity is «f ~ 1/ (ma) . 
Hence, e 2 /(eoVp) ~ 1. This estimate is also valid in 
graphene. Hence, unlike in QED, where chqed = 1/137, 
the coupling constant in graphene is a ~ 1. 

Despite the fact that the coupling constant is of the 
order of unity, a perturbative RG analysis can be applied. 
RG techniques allow us to identify stable fixed points 
of the model, which may be attractive over a broader 
range than the one where a perturbative treatment can 
be rigorously justified. Alternatively, an RG scheme can 
be reformulated as the p rocess of piecew ise integration of 
high energy excitations ( Shankarl . ll994f ). This procedure 
leads to changes in the effective low energy couplings. 
The scheme is valid if the energy of the renormalized 
modes is much larger than the scales of interest. 

The Hartree-Fock correction due to Coulomb interac- 
tions between electrons (given by the diagram in Fig. [37j) 
gives a logarithmic cor rection to the electron self-energy 
([Gonzalez et aLl . ll994h : 



£^(k) = ffcln(^ 



(220) 



where A is a momentum cutoff which sets the range 
of validity of the Dirac equation. This result re- 
mains true even to higher order in perturbation theory 
( Mishchenkol. 12007) and is also obtained in large A^ ex- 
pansions (jRosenstein et all . Il989l . Il99ll : ISonl . 12007?) (N is 
the number of flavors of Dirac fermions), with the only 
modification being the prefactor in (|22Q[) . This result 
implies that the Fermi velocity is renormalized towards 
higher values. As a consequence, the density of states 
near the Dirac energy is reduced, in agreement with the 



general trend of repulsive interactions to induce or in- 
crease gaps. 

This result can be understood from the RG point of 
view by studying the effect of reducing the cut-off from 
A to A — dA and its effect on the effect ive coupling. It 
can b e shown that a obeys the equation ( Gonzalez et all . 
Il994h : 
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Therefore, the Coulomb interaction becomes marginally 
irrelevant. These features are confirmed by a full rela- 
tivistic calculation, although the Fermi velocity cannot , 



obvio usly, surpass the velocity of light l| Gonzalez et all . 
11994 ). This result indicates that stron gly correlated elec - 
tronic phases, such as ferromagnetism ([Peres et aZ.I . l2005h 
and Wigner crystals ( Dahal et all . l2006f ) are suppressed 



in clean graphene. 

A calculation of higher order self-energy terms leads 
to a wavefunction renormalization, and to a finite quasi- 
particle lif etime, which gr o ws linearly with quasiparti- 
cle energy (jGonzalez et all . 1 19941 [1999 ). The wavefunc- 



tion renormalization implies that the quasiparticle weight 
tends to zero as its energy is reduced. A strong cou- 
pling expansion is also possible, assuming that the num- 
ber of electronic flavors justifies an RPA expansion, keep- 



ing o nly electron-hole bubble diagrams (jGonzalez et all . 
1 1999h . This analysis confirms that the Coulomb interac- 
tion is renormalized towards lower values. 

The enhancement in the Fermi velocities leads to a 
widening of the electronic spectrum. This is consistent 
with measurements of the gaps in narrow single wall 
nanotubes, which show deviations from the scaling with 
R , whe re R is the rad i us, ex pected from the Dirac 
equation (jKane and Meld . l2004j) . The linear depen- 



dence of the inverse quasiparticle lifetime with energy is 
consistent with photo-emission experiments in graphite, 



actions ( 


Bostwic 


4 et all. 2007ri Sut 


[awara et all. 120071 


IXu et all 


. 19961: 


Zhou et all. l2006al 


(J). Note that in 



graphite, ba nd structure effects m odify the lifetimes at 



low energies ( Spataru et aLl . I2OO1I ) . The vanishing of the 
quasiparticle peak at low energies can lead to an en- 
ergy dependent renormalization of the interlayer hopping 
( Vozmediano et all . 120021 . |2003[ ). Other thermodynamic 
properties of undoped and doped graphene can also be 
calculated (|Barlas et all . 120071 : IVafekl . 120071 ) . 

Non-perturbative calculations of the effects of the 
long range interactions in undoped graphene show 
that a transition to a gapped phase is also possi- 
ble, when the number of electronic flavors is large 
dKhveshchenkcl l200ll: iKhveshchenko and ShiveR l2006t 
iLuk'vanchuk and Kopelevichl . 120041 ). The broken sym- 
metry phase is similar to the excitonic transition found in 
materials where it becomes favorable to create electron- 
hole pairs that then form bound excitons (excitonic tran- 
sition). 

Undoped graphene cannot have well defined plasmons, 
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as their energies fall within the electron-hole continuum, 
and therefore have a significant Landau damping. At fi- 
nite temperatures, however, thermally excited quasipar- 
ticles screen the Coulomb interactio n, and an ac oustic 
collective charge excitation can exist 

Doped graphene shows a finite density of states at the 
Fermi level, and the long range Coulomb interaction is 
screened. Accordingly, there are collective plasma inter- 
actions near q — > 0, which disperse as u)„ ~ ^/ffi' since 
the sy stem is 2D i Campagnoli and Tosattiri989l : IShungl . 
Il986al fbh . The fact that the electronic states are de- 
scribed by the massless Dirac equation implies that top oc 
n 1 / 4 , where n is the carrier density. The static dielectric 
constant has a continuous de rivativ e at 2fcp unlike in the 
case of the 2D elect r on ga s 1 Andd . 120061% Sarma et all 
120071 : IWunsch et all . l200fih . This fact is associated with 
the suppressed backward scattering in graphene. The 
simplicity of the band structure of graphene allows an- 
alytical calculation of the energy and momentum de- 
pendence of the die lectric function 1 Sarma et all 120071 : 
IWunsch et all l2006h . The screening of the long-range 
Coulomb interaction implies that the low energy quasi- 
particles show a quadratic dependence on ener gy with 
respect to the Fermi energy IjHwang et "aZI. l2007t l. 

One way to probe the strength of the electron- 
electron interactions is via the electronic compress- 
ibility. Measurements of the compressibility using 
a single electron transistor (SET) show very little 
sign of interactions in the system, being well fitted 
by the non- interacting result that, contrary to the 
two-d i mensional electron ga s (2DE G) ([Eisenstein et all 
Il994l : iGiuliani and VignaleL 120051). is positiv ely diver- 
gent ([Martin et all 120071 : IPolini et al\ . l2007h . Bilayer 
graphene, on the other hand, shares characteristics of 
the single layer and the 2DEG with a non-monotonic de- 
pendence of the compre ssibility on the carrier density 
( Kusminskiy et a/1 . l2007l ). In fact, bilayer graphene very 
close to half-filling has been p redicted to be un stable 
towards Wigner crystallization ([Dahal et aZ.I . l2007h . just 
like the 2DEG. Furthermore, according to Hartree-Fock 
calculations, cle an bilayer graphene is unstable towards 
ferromagnetism ( Nilsson et all l2006bh . 



1. Screening in graphene stacks 

The electron-electron interaction leads to the screen- 
ing of external potentials. In a doped stack, the charge 
tends to accumulate near the surfaces, and its distribu- 
tion is determined by the dielectric function of the stack 
in the out-of-plane direction. The same polarizability 
describes the screening of an external field perpendicular 
to the layers, like the one induced by a gat e in electri- 
cally doped systems ( Novoselov et all l2004f l. The self- 
consistent distribution of charge in a biased grap hene bi- 
layer has been studied in ref. I McCannl l2006h . From 
the observed charge distribution and self-consistent cal- 
culations, an estimate of the band structure parameters 



and their rela t ion wi th the induced gap can be obtained 
(|Castro et aZJ . l2007aj) . 

In the absence of interlayer hopping, the polarizability 
of a set of stacks of 2D electron gases can be written as 
a sum of the screening by the individual layers. Using 
the accepted values for the effective masses and carrier 
densities of graphene, this scheme gives a first approx- 
imati on to screening in graphite ( Visscher and Falicovl . 
1 1 9 7 lh . The screening length in t he out of plane direc - 
tion is of about 2 graphene layers ( Morozov et all [20051 ) . 
This model is easily generalizable to a s tack of semimet- 
als d escribed by the 2D Dirac equation ( Gonzalez et all 
l200lh . At half filling, the screening length in all direc- 
tions diverges, and the screening effects are weak. 

Interlayer hopping modifies t his picture sig nificantly. 
The hopping leads to coherence ( Guineal . l2007t i. The out 
of plane electronic dispersion is similar to that of a one 
dimensional conductor. The out of plane polarizability 
of a multilayer contains intra- and interband contribu- 
tions. The subbands in a system with the Bernal stack- 
ing have a parabolic dispersion, when only the nearest 
neighbor hopping terms are included. This band struc- 
ture leads to an interband susceptibility described by a 
sum of terms like those in l|228|) , which diverges at half- 
filling. In an infinite system, this divergence is more pro- 
nounced at k± = 7r/c, that is, for a wave vector equal 
to twice the distance between layers. This effect greatly 
enhances Friedel like oscillations in the charge distribu- 
tion in the out of plane direction, which can lead to the 
changes i n the sign of the charge in neighboring layers 
([Guinea! . l2007t ). Away from half-filling a graphene bi- 
layer behaves, from the point of view of screening, in a 
way ve ry similar to the 2DEG ( Wang and Chakrabortyl 
l2007bh . 



C. Short range interactions 

In this section we discuss the effect of short range 
Coulomb interactions on the physics of graphene. The 
simplest carbon system with a hexagonal shape is the 
benzene molecule. The value of the Hubbard interaction 
among 7r-electrons was, for this sys tem, computed long 
ago by Parr et al. ( Parr et aZl . ll950h . yielding a value of 
U = 16.93 eV. For comparison purposes, in polyacetylene 
the value for the Hubbard inter action is U ~10 eV an d 
the hopping energy is t f=a 2.5 eV ( Baeriswvl et «/.Ul986f l. 
These two examples just show that the value of the on- 
site Coulomb interaction is fairly large for 7r— electrons. 
As a first guess for graphene, one can take U to be of the 
same order as for polyacethylene, with the hopping inte- 
gral t ~ 2.8 eV. Of course in pure graphene the electron- 
electron interaction is not screened, since the density of 
states is zero at the Dirac point, and one should work 
out the effect of Coulomb interactions by considering the 
bare Coulomb potential. On the other hand, as we have 
seen before, defects induce a finite density of states at the 
Dirac point, which could lead to an effective screening of 
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the long-range Coulomb interaction. Let us assume that 
the bare Coulomb interaction is screened in graphene and 
that Coulomb interactions are represented by the Hub- 
bard interaction. This means that we must add to the 
Hamiltonian ((5j) a term of the form: 



Hu = U 
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The simplest question one can ask is whether this system 
shows a tendency toward some kind of magnetic order 
driven by the interacti on U. Within the s implest Hartree- 
Fock approximation ( Peres et~al\ . l2004h . the instability 
line toward ferromagnetism is given by: 



U F (n) 
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which is nothing but the Stoner criterion. Simila r results 
are o btained in more sophisticated calculations 1 Herbutl . 
l2006h . Clearly, at half-filling the value for the density of 
states is p(0) = and the critical value for Uf is arbi- 
trarily large. Therefore we do not expect a ferromagnetic 
ground state at the neutrality point of one electron per 
carbon atom. For other electronic densities, p(p) be- 
comes finite producing a finite value for Uf- We note 
that the inclusion of *' does not change these findings, 
since the density of states remains zero at the neutrality 
point. 

The lin e toward an antife rromagnetic ground state is 
given by (jPeres et all . l2004h 



Uaf(p) = 
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where E + (k) is given in (jgj. T his result give s a fi- 
nite V af at the neutrali ty point ( Martelo et all . Il997t 
ISorella and Tosattil . fl992T ) : 



U AF {0) = 2.23*. 



(225) 



Quantum Monte Carlo c alculations i Paiva et all . l2005t 
ISorella and Tosattil . 11992^ 1 . raise however its value to: 



U AF (0)~5t, 



(226) 



Taking for graphene the same value for U as in polyacety- 
lene and t = 2.8 eV, one obtains U/t ~ 3.6, which put the 
system far from the transition toward an antiferromagnet 
ground state. Yet another possibility is th at the system 
may be in a sort of a quantum spin liqui d llLee and Lee , 
l2005h (as originally proposed by Pauling ijPaulingl . Il972 ) 
in 1956) since mean field calculations give a critical value 
for U to be of the order of U/t ~ 1.7. Whether this type 
of ground state really exists and whether quantum fluc- 
tuations pushes this value of U toward larger values is 
not known. 



1. Bilayer graphene: exchange 

The exchange interaction can be large in an unbiased 
graphene bilayer with a small concentration of carriers. It 
was shown previously that the exchange contribution to 
the electronic energy of a single grap hene layer does no t 
lead to a ferromagnetic instability 1 Peres et all . 120051 ). 
The reason for this is a significant contribution from the 
interband exchange, which is a term usually neglected 
in doped semiconductors. This contribution depends on 
the overlap of the conduction and valence wavefunctions, 
and it is modified in a bilayer The interband exc hange 
energy is reduced in a bilayer ( Nilsson et all . l2006d) . and 
a positive contribution that depends logarithmically on 
the bandwidth in graphene is absent in its bilayer. As a 
result, the exchange energy becomes negative, and scales 
as n 3 / 2 , where n is the carrier density, similar to the 
2DEG. The quadratic dispersion at low energies implies 
that the kinetic energy scales as n 2 , again as in the 2DEG. 
This expansion leads to: 
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Writing rij = (n + s)/2,nj = (n — s)/2, where s is the 
magnetization, (|227p predicts a second order transition to 
a ferromagnetic state for n = (4e 4 * 2 )/(8l7r 3 ii|,en). Higher 
order corrections to l|227p l ead to a first order tr ansition 
at slightly higher densities ijNilsson et all . l2006d) . For a 
ratio 71/70 ~ 0.1, this analysis implies that a graphene 
bilayer should be ferromagnetic for carrier densities such 
that \n\ < 4 x 10 10 crrr 2 . 

A bilayer is also the unit cell of Bernal graphite, and 
the exchange instability can also be studied in an infi- 
nite system. Taking into account nearest neighbor in- 
terlayer hopping only, bulk graphite should also show 
an exchange instability at low doping. In fact, there is 
some experimental evidence for a ferr omagnetic insta- 
bility in strongly disordered graphite llEsquinazi et all . 
12002 . l2003t iKopelevich and Esquinazl 12000 1. 

The analysis described above can be extended to the 
biased bilayer, where a gap separate s the conduction and 
valence bands ( Stauber et a/.l . l200^ 1. The analysis of this 
case is somewhat different, as the Fermi surface at low 
doping is a ring, and the exchange interaction can change 
its bounds. The presence of a gap reduces further the 
mixing of the valence and conduction band, leading to an 
enhancement of the exchange instability. At all doping 
levels, where the Fermi surface is ring shaped, the biased 
bilayer is unstable towards ferromagnetism. 



2. Bilayer graphene: short range interactions 

The band structure of a graphene bilayer, at half fill- 
ing, leads to logarithmic divergences in different response 
functions at q = 0. The two parabolic bands that are 
tangent at k = lead to a susceptibility which is propor- 
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Figure 38 (Color online) Sketch of the expected magnetiza- 
tion of a graphene bilayer at half-filling. 



tional to: 



d 2 k 



1 



oc log 



A 



* lm<A ' - (4A)|k| 2 " ~° V VMM 

(228) 

where A ~ \A 2 / u f i s a high momentum cutoff. These 
logarithmic divergences are similar to the ones which 
show up when the Fermi surface of a 2 D metal is near a 
saddl e point in the dispersion relation ( Gonzalez et ali . 
Il996h . A full treatment of th ese divergences requires a 



RG approach ( Shankarl . Il994h . Within a simpler mean 
field treatment, however, it is easy to notice that the di- 
vergence of the bilayer susceptibility gives rise to an in- 
stability towards an antiferromagnetic phase, where the 
carbon atoms which are not connected to the neighbor- 
ing layers acquire a finite magnetization, while the mag- 
netization of the atoms with neighbors in the contiguous 
layers remain zero. A scheme of the expected ordered 
state is shown in Fig. [38l 



D. Interactions in high magnetic fields 

The formation of Landau levels enhances the effect of 
interactions due to the quenching of the kinetic energy. 
This effect is most pronounced at low fillings, when only 
the lowest levels are occupied. New phases may appear 
at low temperatures. We consider here phases differ- 
ent from the fractional quantum Hall effect, which has 
not been observed in graphene so far. The existence of 
new phases can be inferred from the splitting of the val- 
ley or spin degeneracy of the Landau level s, which can 
be obs e rved in spectros c opy m easurements ( Jiang et~al\ . 
l2007al : ISadowski et ali . l2006h , o r in the app earance 



2UUYat ISadowski et all I2UUO ). o r m trie app 
of new quantum Hall plateaus (lAbanin et ali. 



Gicsbc 



2007bl : 



Zhang et al. 



s et ali. 2007 : Goswami et ali . 2007 : Jiang et 



20061 ) 



2007b; 



Interactions can lead to new phases when their effect 
overcomes that of disorder. An analysis of the com- 
petit ion between disorder and inter actions is found in 
ref. ( Nomura and MacDonaldl . l2007h . The energy split- 
ting of the different broken symmetry phases, in a clean 
system, is determined by lattice effects, so that it is re- 
duced by factors of order a/Is, where a is a length of 
the ord er of the lattice sp a cing, a nd Ir i s the magnetic 
length dAlicea and Fished . 120061 . 120071 : iGoerbig et~ali . 
120061 : 1 Wang et ali . l2007h . The combination of disorder 



and a magnetic field may also lift the degeneracy be- 
tween the two valley s, favoring valley polarized phases 
(| Abanin et afl l20Q7ah . 

In addition to phases with enhanced ferromag- 
netism or with broken valley symmetry, interactions 
at high magnetic fie l ds ca n lead to excitonic instabil- 

and Wigner crystal phases 
When only the n = state 



i ties ( Gusvni n et al. 
(jZhang and Jogleka 



2006) 



2007) 



is occupied, the Landau levels have all their weight in a 
given sublattice. Then, the breaking of valley degener- 
acy can be as sociated with a charge de nsity wave, which 
opens a gap ( Fuchs and Lederen . 120071 ). It is interesting 
to note t hat in these phases new collective e xcitations are 
possible (jDoretto and Morais Smit hi . [20071 ). 

Interactions modify the edge states in the quantum 
Hall regime. A novel phase can appear when the n = is 
the last filled level. The Zeeman splitting shifts the elec- 
tron and hole like chiral states, which disperse in opposite 
directions near the boundary of the sample. The result- 
ing level crossing between an electron like level with spin 
anti-parallel to the field, and a hole like level with spin 
parallel to the fiel d, may lead to Luttinger liquid features 
in the edge states ( Abanin et a/J . [2007bHFertig and BreyL 
l2006h . 



VI. CONCLUSIONS 

Graphene is a unique system in many ways. It is truly 
2D, has unusual electronic excitations described in terms 
of Dirac fermions that move in a curved space, it is an in- 
teresting mix of a semiconductor (zero density of states) 
and a metal (gaplessness) , and has properties of soft mat- 
ter. The electrons in graphene seem to be almost insen- 
sitive to disorder and electron-electron interactions and 
have very long mean free paths. Hence, graphene's prop- 
erties are rather different from what is found in usual 
metals and semiconductors. Graphene has also a very 
robust but flexible structure with unusual phonon modes 
that do not exist in ordinary 3D solids. In some sense, 
graphene brings together issues in quantum gravity and 
particle physics, and also from soft and hard condensed 
matter. Interestingly enough, these properties can be 
easily modified with the application of electric and mag- 
netic fields, addition of layers, by control of its geometry, 
and chemical doping. Moreover, graphene can be directly 
and relatively easily probed by various scanning probe 
techniques from mesoscopic down to atomic scales, be- 
cause it is not buried inside a 3D structure. This makes 
graphene one of the most versatile systems in condensed 
matter research. 

Besides the unusual basic properties, graphene has 
the potential for a large number of applications 

fro m ch emical sensors 
2007|l to transistors 



Geim and Novoseloyl. 1 2007) 



Chen et ali. 2007c : Schedin et al 



Nilsson et ali . l2007bHOostinga et all l2007h . Graphene 



can be chemically and/or structurally modified in order 
to change its functionality and henceforth its potential 
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applications. Moreover, graphene can be easily obtained 
from graphite, a material that is abundant on earth's 
surface. This particular characteristic makes graphene 
one of the most readily available materials for basic re- 
search since it frees economically challenged research in- 
stitutions in developing countries from the dependence of 
expensive sample growing techniques. 

Many of graphene's properties are currently subject of 
intense research and debate. The understanding of the 
nature of the disorder and how it affects the transport 
properties (a problem of fundamental importance for ap- 
plications), the effect of phonons on electronic transport, 
the nature of electron-electron interactions and how they 
modify its physical properties are research areas that are 
still in their infancy. In this review, we have touched 
only the surface of a very deep sea that still has to be 
explored. 

Whereas hundreds of papers have been written on 
monolayer graphene in the last few years, only a small 
fraction actually deals with multilayers. The majority of 
the theoretical and experimental efforts have been con- 
centrated on the single layer, perhaps because of its sim- 
plicity and the natural attraction that a one atom thick 
material, which can be produced by simple methods in 
almost any laboratory in the world, creates for human 
imagination. Nevertheless, few layer graphene is equally 
interesting and unusual with a technological potential 
perhaps bigger than the single layer. Indeed, the the- 
oretical understanding and experimental exploration of 
multilayers is far behind the single layer. This is a fertile 
and open field of research for the future. 

Finally, we have focused entirely on pure carbon 
graphene where the band structure is dominated by the 
Dirac description. Nevertheless, chemical modification 
of graphene can lead to entirely new physics. Depend- 
ing on the nature of chemical dopants and how they are 
introduced into the graphene lattice (adsorption, substi- 
tution, or intercalation) the results can be many. Small 
concentrations of adsorbed alkali metal can be used to 
change the chemical potential while adsorbed transition 
elements can lead to strong hybridization effects that af- 
fect the electronic structure. In fact, the introducion of d- 
and f-electron atoms in the graphene lattice may produce 
a significant enhancement of the electron-electron inter- 
actions. Hence, it is easy to envision a plethora of many- 
body effects that can be induced by doping and have 
to be studied in the context of Dirac electrons: Kondo 
effect, ferromagnetism, antiferromagnetism, charge and 
spin density waves. The study of chemically induced 
many-body effects in graphene would add a new chap- 
ter to the short but fascinating history of this material. 
Only future will tell but the potential for more amaze- 
ment is lurking on the horizon. 
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